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ABSTRACT 

Planck has produced detailed all-sky observations over nine frequency bands between 30 and 857 GHz. These observations allow robust recon- 
struction of the primordial cosmic microwave background (CMB) temperature fluctuations over nearly the full sky, as well as new constraints on 
Galactic foregrounds, including thermal dust and line emission from molecular carbon monoxide (CO). This paper describes the component sepa- 
ration framework adopted by Planck for many cosmological analyses, including CMB power spectrum determination and likelihood construction 
on large angular scales, studies of primordial non-Gaussianity and statistical isotropy, the integrated Sachs-Wolfe effect (ISW), and gravitational 
lensing, and searches for topological defects. We test four foreground-cleaned CMB maps derived using qualitatively different component separa- 
tion algorithms. The quality of our reconstructions is evaluated through detailed simulations and internal comparisons, and shown through various 
tests to be internally consistent and robust for CMB power spectrum and cosmological parameter estimation up to ^ = 2000. The parameter 
constraints on ACDM cosmologies derived from these maps are consistent with those presented in the cross-spectrum based Planck likelihood 
analysis. We choose two of the CMB maps for specific scientific goals. We also present maps and frequency spectra of the Galactic low-frequency, 
CO, and thermal dust emission. The component maps are found to provide a faithful representation of the sky, as evaluated by simulations, with 
the largest bias seen in the CO component at 3 %. Eor the low-frequency component, the spectral index varies widely over the sky, ranging from 
about p - -Aio -2. Considering both morphology and prior knowledge of the low frequency components, the index map allows us to associate a 
steep spectral index {p < -3.2) with strong anomalous microwave emission, corresponding to a spinning dust spectrum peaking below 20 GHz, a 
flat index ofp > -2.3 with strong free-free emission, and intermediate values with synchrotron emission. 

Key words, cosmology: observations 
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1. Introduction 

This paper, one of a set associated with the 2013 release of 
data from the Planck^ mission (Planck Collaboration I 2013), 
describes the component separation techniques applied to the 
Planck data to produce maps of the cosmic microwave back- 
ground (CMB) temperature anisotropics (see Fig. 1) and of dif- 
fuse foregrounds. 

The sky at millimetre and sub-millimetre wavelengths con- 
tains a wealth of cosmological and astrophysical information. 
Accessing it is an inversion process, known as component sepa- 
ration, to extract the sources of emission contributing to a set of 
maps observed at different frequencies. Planck gives us a power- 
ful data set to unlock new information in this manner by observ- 
ing the entire sky from 30 to 857 GHz in nine frequency bands 
at higher angular resolution and sensitivity than its predecessors. 
Accurate and detailed component separation is a central objec- 
tive of the mission. 

We divide the foregrounds into two distinct categories: 
diffuse emission from the Galaxy and compact sources. The 
Galactic foregrounds are the principal source of contamina- 
tion of the CMB on large angular scales, with fluctuation 
power decreasing roughly as a power law towards higher mul- 
tipoles (Bennett et al. 2003). They are dominated by syn- 
chrotron, free-free and Anomalous Microwave Emission (AME, 
ascribed to spinning dust grains, at frequencies below 70 GHz) 
and by rotational line emission from carbon monoxide (CO) 
molecules and thermal dust emission at frequencies above 
100 GHz. Extragalactic foregrounds, on the other hand, domi- 
nate the small-scale contamination of the CMB. They arise from 
discrete, individually detectable compact sources and the collec- 
tive emission from unresolved radio and infrared (IR) sources, 
and also from the Sunyaev-Zeldovich (SZ) effect in galaxy clus- 
ters (Planck Collaboration XVIII 2011; Planck Collaboration 
XXVIII 2013; Planck Collaboration XXIX 2013). 

In the Planck analyses, these foregrounds are dealt with in 
a variety of ways. At the power spectrum and likelihood level, 
the extragalactic foregrounds are modelled with parameterized 
power spectra, appropriate to their statistical isotropy, over re- 
gions restricted to low Galactic emission (Planck Collaboration 
XV 2013). Component separation as described in the present pa- 
per aims at removing Galactic emission to produce CMB maps 
covering the largest possible sky area for studies of the large- 
scale properties and higher-order statistics of the CMB. In ad- 
dition, this component separation provides a reconstruction of 
the diffuse emission from our Galaxy. Detailed studies of spe- 
cific extragalactic foregrounds, such as the cosmic infrared back- 
ground (CIB) (Planck Collaboration XVIII 2013) and the dif- 
fuse Sunyaev-Zeldovich (SZ) signal (Planck Collaboration XXI 
2013), employ methods tailored to their particular needs. 

Building on previous work (Leach et al. 2008), we approach 
CMB extraction with a philosophy designed to ensure robust- 
ness by applying four distinct algorithms based on two different 
methodologies. The first avoids any assumptions concerning the 
foregrounds and relies solely on a minimum variance criterion 
for the data component possessing a blackbody spectrum (i.e., 
the CMB), while the second methodology relies on parametric 
modelling of the foregrounds in either real or harmonic space. 

^ Planck (http://www.esa.int/Planck) is a project of the 
European Space Agency (ESA) with instruments provided by two sci- 
entific consortia funded by ESA member states (in particular the lead 
countries France and Italy), with contributions from NASA (USA) and 
telescope reflectors provided by a collaboration between ESA and a sci- 
entific consortium led and funded by Denmark. 



We evaluate the performance of these component separation al- 
gorithms through detailed simulations, and we examine the ro- 
bustness of the recovered CMB maps by comparing them, their 
power spectra, and their resulting cosmological constraints. As a 
diagnostic, we also briefly examine their higher-order statistics. 

The CMB results presented in this work serve a number ap- 
plications. We use the real-space modelling to produce a clean 
CMB map and power spectra on large angular scales, where 
diffuse Galactic emission is the main contaminant, to construct 
the likelihood function at low multipoles; this is then combined 
with the high multipole likelihood function that models extra- 
galactic foregrounds with power spectra (Planck Collaboration 
XV 2013). The high resolution CMB maps are used as a check 
on primary cosmological constraints (see below), for lensing 
studies (Planck Collaboration XVII 2013), studies of the inte- 
grated Sachs-Wolfe effect (Planck Collaboration XIX 2013), of 
the isotropy of the CMB (Planck Collaboration XXIII 2013), 
of non-Gaussian statistics (Planck Collaboration XXIV 2013), 
in searches for topological defects (Planck Collaboration XXV 
2013), and for examination of the geometry and topology of the 
Universe (Planck Collaboration XXVI 2013). 

In addition, we present maps of diffuse Galactic emission 
divided into low- and high-frequency components, as well as a 
molecular CO component. We judge the adequacy of this recon- 
struction through simulations and by comparison with known 
properties of the diffuse Galactic foregrounds. 

The paper is organized as follows. In Sect. 2 we discuss the 
expected sources of sky emission over the Planck frequency in- 
terval and how they are modelled. Then in Sect. 3 we detail the 
overall approach and introduce the four component separation 
methods. In Sect. 4 we present the Planck data set and pre- 
processing procedure, and we describe our simulations. This is 
followed by a presentation of the derived CMB maps and their 
characterization in Sect. 5. Section 6 is dedicated to power spec- 
tra and cosmological parameter constraints obtained from these 
maps, and Sect. 7 to studies of higher-order statistics. Section 8 
presents a reconstruction of the difl'use Galactic foregrounds, and 
Sect. 9 concludes. We relegate details of the algorithms to appen- 
dices. 

2. The sky at Planck frequencies 

The properties of Galactic emission vary significantly across 
the Planck frequency range from 30 to 857 GHz. At frequen- 
cies below 70 GHz, the dominant radiation processes are: syn- 
chrotron emission from cosmic ray electrons interacting with 
the Galactic magnetic field (e.g., Haslam et al. 1982; Reich & 
Reich 1988; Broadbent et al. 1989; Davies et al. 1996; Platania 
et al. 2003; Bennett et al. 2003; Gold et al. 2011); thermal 
Bremsstrahlung (or free-free emission) from electron-electron 
and electron-ion scattering (e.g., Banday et al. 2003; Dickinson 
et al. 2003; Davies et al. 2006; Ghosh et al. 2012; Planck 
CoUaboration Int. XII 2013; Planck Coflaboration XX 2011); 
and AME from dust grains (Kogut 1996; Leitch et al. 1997; 
Banday et al. 2003; Lagache 2003; de Oliveira-Costa et al. 2004; 
Finkbeiner et al. 2004; Davies et al. 2006; Bonaldi et al. 2007; 
Dobler & Finkbeiner 2008; Miville-Deschenes et al. 2008; Ysard 
et al. 2010; Gold et al. 2011; Planck Collaboration XX 2011), 
possibly due to their rotational line emission (Draine & Lazarian 
1998; AH-Haimoud et al. 2009; Ysard & Verstraete 2010; Hoang 
& Lazarian 2012). Over the frequency range covered by Planck, 
both synchrotron and free-free spectra are well approximated by 
power laws in brightness temperature, oc y^, with the syn- 
chrotron index, y^synch, ranging from -3.2 to -2.8 (Davies et al. 
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Fig. 1: Foreground-cleaned CMB maps derived by Commander-Ruler, NILC, SEVEM and SMICA. Note that the SMICA map has been 
filled in smoothly inside a 3 % Galactic mask. 



1996) and the free-free index, ySff, lying between -2.2 and -2.1. 
Less is known about the AME spectrum, but spinning dust mod- 
els with a spectrum peaking at frequencies below 20 GHz (in 
brightness temperature units) adequately describe current obser- 
vations^. Above the peak, the spectrum appears consistent with 
a power-law (e.g., Banday et al. 2003; Davies et al. 2006; Dobler 
& Finkbeiner 2008; Ghosh et al. 2012). In addition to these three, 
the existence of a fourth low-frequency foreground component, 
known as the "Galactic haze", has been claimed, possibly due to 
a hard- spectrum synchrotron population near the Galactic cen- 
tre (e.g., Finkbeiner 2004; Dobler & Finkbeiner 2008; Pietrobon 
et al. 2012; Planck Collaboration Int. IX 2013). 

At frequencies higher than 100 GHz, thermal dust emission 
dominates over most of the sky and is commonly described 
by a modified blackbody spectrum with power-law emissivity, 
6y oc y^d^ and temperature, T^. Both the temperature and spectral 
index, vary spatially. Prior to Planck, the best-fitting single 
component dust model had a temperature ^ 18 K and spec- 
tral index ^ 1.7 (Finkbeiner et al. 1999; Bennett et al. 2003; 
Gold et al. 2011), although there is evidence of flattening of the 
spectral index from around 1.8 in the far-infrared to 1.55 in the 



^ Note that we adopt brightness temperature for AME in this pa- 
per, while many other publications adopt flux density. When compar- 
ing peak frequencies, it is useful to note that that a spectrum that has 
a maximum at 30 GHz in flux density peaks at 17 GHz in brightness 
temperature. 



microwave region (Planck Collaboration 2012), the interpreta- 
tion of which is still under study. 

In addition to these difl'use Galactic components, extra- 
galactic emission contributes at Planck frequencies. In partic- 
ular, a large number of radio and far-infrared (FIR; Planck 
Collaboration XIII 2011) galaxies, clusters of galaxies and the 
Cosmic Infrared Background (CIB; Planck Collaboration XVIII 
2011) produce a statistically isotropic foreground, with fre- 
quency spectra well approximated by models similar to those ap- 
plicable to the Galactic foregrounds (modified blackbody spec- 
tra, power laws, etc.). Except for a frequency-dependent absolute 
off'set, which may be removed as part of the overall off'set re- 
moval procedure, these extragalactic components are therefore 
typically absorbed by either the low-frequency or thermal dust 
components during component separation. No special treatment 
is given here to extragalactic foregrounds, beyond the masking 
of bright objects. Dedicated scientific analyses of these sources 
are described in detail in Planck Collaboration XVIII (2011), 
Planck Collaboration XXVIII (2013), and Planck Collaboration 
XXIX (2013). In the Planck likelihood, extragalactic sources are 
modelled in terms of power spectrum templates at high I (Planck 
Collaboration XV 2013). 

Other relevant sources include emission from molecular 
clouds, supernova remnants, and compact H ii regions inside our 
own Galaxy, as well as the thermal and kinetic SZ eff'ects, due to 
inverse Compton scattering of CMB photons off' free electrons 
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Table 1: Overview and comparison of component separation algorithms. 



Characteristic 


C oimnand e r - Ru 1 e r 


NILC 


SEVEM 


SMICA 


Method 


Bayesian parameter 


Internal linear 


Internal template 


Spectral parameter 




estimation 


combination 


fitting 


estimation 


Domain 


Pixel 


Needlet 


Pixel 


Spherical harmonic 


Channels [GHz] 


30-353 


44-857 


30-857 


30-857 


Effective beam FWHM [arcmin] 


-7.4 


5.0 


5.0 


5.0 


/ 


none 


3200 


3100 


4000 






Fig. 2: Combined Galactic (CG) emission masks for the Planck 
data, corresponding to sky fractions of 20, 40, 60, 70, 75, 80, 90, 
97, and 99 %. The masks are named CG20, etc. 



in ionized media. Planck provides new and important informa- 
tion on all these processes, as described both in the following 
and in the companion papers Planck Collaboration XIII (2013), 
Planck Collaboration XXI (2013), and Planck Collaboration 
(2013b). In particular, Planck'^ frequency range, angular reso- 
lution and sensitivity make it a powerful probe of thermal dust, 
resulting in new and tight constraints on dust temperature and 
emissivity. The same frequencies also allow extraction of the 
first ever full-sky maps of the emission resulting from the CO 
/=1^0, /=2^1 and /=3^2 rotational transitions at 115, 230 
and 345 GHz, respectively (Planck Collaboration XIII 2013). 

The focus of this paper is to reconstruct the CMB 
anisotropics over a large sky fraction, exploiting only the Planck 
frequency bands. We also present a detailed reconstruction of the 
thermal dust emission at high frequencies, as well as CO emis- 
sion lines. At low frequencies and over the region used for CMB 
analysis, the total foreground contribution is well approximated 
by a single power law (see Sect. 8). We therefore model the sum 
of all low-frequency foregrounds by a power law with spatially 
varying spectral index whose numerical value in any pixel re- 
sults from the influence of the dominant foreground component 
at that location. The full analysis of diflTuse foregrounds, using 
ancillary data to resolve the individual components at low fre- 
quencies, will be presented in a forthcoming publication. 



3. Approach to component separation 

The rich content of the Planck data encourages application of 
several component separation techniques. We consider four, as 
summarized in Table 1, which we classify according to one 
of two diflTerent general methodologies. The first makes min- 
imal assumptions concerning the foregrounds and seeks only 
to minimize the variance of the CMB, i.e., the sky component 
possessing a blackbody spectrum. We implement this approach 



with a needlet (wavelet on the sphere) version of the internal 
linear combination (ILC) algorithm (NILC; Delabrouille et al. 
2009), and also with a template-based method to remove fore- 
ground contamination from the CMB -dominant bands. These 
foreground templates are constructed from the lowest and high- 
est frequency channels (Spectral Estimation Via Expectation 
Maximization, SEVEM; Femandez-Cobos et al. 2012). 

The second methodology uses parametric modelling of the 
foregrounds. In our real space implementation, we explore 
model parameters through Bayesian parameter estimation tech- 
niques, fitting a parametric signal model per pixel (Coimnander; 
Eriksen et al. 2006, 2008); a similar implementation is pre- 
sented by Stompor et al. (2009). To estimate spectral indices 
robustly in pixel space, this procedure requires identical angu- 
lar resolution across all frequencies included in the analysis, and 
is therefore limited in resolution by the 30 GHz LEI channel. 
However, this is sufliicient to generate the low-resolution CMB 
map and power spectrum samples required for the low multipole 
part of the Planck likelihood function for cosmological param- 
eters (Planck Collaboration XV 2013). To produce full resolu- 
tion maps, we use the resulting low-resolution spectral parameter 
samples to solve for the component amplitudes, in an extension 
to the method known as Ruler (we refer to the combined method 
as Commander-Ruler, or C-R). In our fourth technique, we im- 
plement a CMB -oriented parametric approach that fits the am- 
plitude and spectral parameters of CMB and foregrounds in the 
harmonic domain (Spectral Matching Independent Component 
Analysis, SMICA; Cardoso et al. 2008). 

Details of each algorithm are given in the appendices. We 
now turn to their application to the data and evaluate their per- 
formance using simulations. 



4. Data, simulations and masl<s 

We use the data set from the first 15.5 months of Planck obser- 
vations, corresponding to 2.6 sky surveys, from both the Low 
Erequency Instrument (LEI) and High Erequency Instrument 
(HEI). The primary inputs for component separation are the fre- 
quency channel maps, including half-ring maps, bandpasses, and 
beam characteristics; a full description of these products is given 
in Planck Collaboration II (2013) and Planck Collaboration VI 
(2013). No special corrections are made for zodiacal light emis- 
sion (ZLE; Planck Collaboration VI 2013) in the analyses pre- 
sented here. The ZLE is not stationary on the sky, since it de- 
pends on Planck'^ position and scanning strategy. Therefore the 
frequency maps contain a projected version of the emission av- 
eraged over the nominal mission. Despite this, a series of ex- 
ploratory analyses showed that our algorithms naturally correct 
for this component within their existing model space. It was also 
found that larger CMB residuals were induced when applying a 
correction based on a ZLE model than when applying no correc- 
tion, most likely due to uncertainties in the model itself. 
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Fig. 3: Summary of Component Separation (CS) confidence 
masks. Each pixel is encoded in terms of a sum in which 
Coimnander-Ruler equals 1 (light blue), NILC equals 2 (dark 
red), SEVEM equals 4 (yellow), and SMICA equals 8 (fight red). 
The masks are named CS-CR75, CS-NILC93, CS-SEVEM76, 
and CS-SMICA89, respectively, reflecting their accepted sky 
fraction. The union mask (U73), used for evaluation purposes 
in this paper, removes all coloured pixels. 



To evaluate and validate our algorithms, we analyse a large 
suite of realistic simulations, the so-called Full Focal Plane 
(FFP) simulations, based on detailed models of the instrument 
and sky. The version used for this data release is denoted FFP6, 
and is described in Planck Collaboration ES (2013). The simula- 
tion procedure generates time streams for each detector, incorpo- 
rating the satellite pointing, the individual detector beams, band- 
passes, noise properties, and data flags, and then produces sim- 
ulated frequency channel maps through the mapmaking process. 
For the input sky, we use the Planck Sky Model (PSM), which 
includes the CMB, difl'use Galactic emission (synchrotron, free- 
free, thermal dust, AME, and molecular CO lines), and com- 
pact sources (thermal and kinetic SZ eff'ects, radio sources, in- 
frared sources, the CIB, and ultra-compact H ii regions). The pre- 
launch version of the PSM is described by Delabrouille et al. 
(2012), and has been modified for the present work as described 
in Planck Collaboration ES (2013). Each FFP data set consists of 
three parts: the simulated observations, Monte Carlo realizations 
of the CMB, and Monte Carlo realizations of the instrumental 
noise. 

For both the data and the simulations, we reconstruct the 
CMB and foregrounds from the full frequency channel maps and 
the corresponding half-ring maps, which are made from the data 
in the first half or second half of each stable pointing period. The 
half-ring maps can be used to obtain an estimate of the noise in 
each channel by taking half of the diff'erence between the two 
maps, thereby normalizing the noise level to that of the full map. 
This is referred to as the half-ring half-diff'erence (HRHD) map. 
The signals fixed to the sky will be cancelled leaving only the 
noise contribution. The HRHD map can be treated as a realiza- 
tion of the same underlying noise processes and it can be used to 
estimate the power spectrum, and other properties, of the noise. 
If there are noise correlations between the half-ring maps, then 
the estimates of the noise properties thus obtained can be biased. 
This is the case for HFI channels; the cosmic ray glitch removal 
(Planck Coflaboration VI 2013; Planck Coflaboration X 2013) 
induces correlations that lead to the noise power spectrum being 
underestimated by a few percent at high £ when using the HRHD 
maps. 




Fig. 4: Beam transfer functions of the four foreground-cleaned 
CMB maps. 




/iK 10 



Fig. 5: Standard deviation between the four foreground-cleaned 
CMB maps. All maps have been downgraded to a HEALPix res- 
olution of /"/side = 128. The diff'erences are typically less than 
5 yuK at high Galactic latitudes, demonstrating that the maps are 
consistent over a large part of the sky. 



Prior to processing the data through each component sep- 
aration pipeline, we define masks for the point sources and 
bright Galactic regions. Point source masking is based on the 
source catalogues obtained by filtering the input sky maps 
with the Mexican Hat Wavelet 2 (MHW2) filter and apply- 
ing a 4(T threshold for the LEI bands and a 5(T threshold 
for the HFI bands (Planck Collaboration XVIII 2011; Planck 
Collaboration XXVIII 2013). The mask radius of each source 
is diff'erent for the LFI and HFI. Due to the large beam size 
of LFI channels, we define a variable masking radius for each 
source according to its signal-to-noise ratio (S/N) as r = 
(2 log(A/m))^/^/(2 ^2 log 2)xFWHM, where r is the radius, A is 
the S/N, and m is the maximum amplitude (given in units of the 
background noise level) allowed for the tail of unmasked point 
sources; we set m = 0.1, which is a compromise between mask- 
ing the source tails and minimizing the number of masked pixels. 
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C-R - NILC C-R - SEVEM 






C-R - SMICA 



NILC - SEVEM 





NILC - SMICA 



SEVEM - SMICA 






-30 



30 



Fig. 6: Pairwise differences between foreground-cleaned CMB maps. All maps have been downgraded to a HEALPix resolution of 
/"/side = 128 to show the large-scale differences. The line-like discontinuities in the differences involving SEVEM is due to the two 
different regions used in this algorithm to clean the sky (see Appendix C for details). 



For HFI, the mask radius around each source is 1.27xFWHM, 
using the average FWHM obtained from the effective beams. 



respond to sky fractions of 20, 40, 60, 70, 75, 80, 90, 97, and 
99 %, and are named CG20, etc. 



A basic set of Galactic masks is defined as follows. We sub- 
tract a CMB estimate from the 30 and 353 GHz maps, mask 
point sources, and smooth the resulting maps by a Gaussian with 
FWHM of 5°. We then threshold and combine them, generating 
a series of masks with different amounts of available sky. The 
resulting combined Galactic (CG) masks, shown in Fig. 2, cor- 



5. CMB Maps 

We begin the discussion of our results by presenting the 
foreground-cleaned CMB maps. These maps are shown in Fig. 1 
for each of the four component separation algorithms. Already 
from this figure it is clear that the wide frequency coverage and 
high angular resolution of Planck allow a faithful reconstruction 
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Fig. 7: CMB residual maps from the FFP6 simulation. A monopole determined at high Galactic latitude has been subtracted from 
the maps, and they have been downgraded to a HEALPix resolution of A^side = 128 to show the large-scale features. The residuals 
presented here provide a conservative estimate of those expected in the data (see text for details). 



of the CMB field over most of the sky. The fluctuations appear 
visually consistent with the theoretical expectation of a Gaussian 
and isotropic signal everywhere except inside a small band very 
close to the Galactic plane^ 

Each CMB map is accompanied by its own confidence mask 
outside which the corresponding solution is considered statis- 
tically robust, shown in Fig. 3; for a definition of each mask, 
see Appendices A-D. Accepted sky fractions are 75, 93, 76, and 
89%, respectively, for Coimnander-Ruler, NILC, SEVEM, and 
SMICA. These masks are denoted CS-CR75, CS-NILC93, CS- 
SEVEM76, CS-SMICA89, respectively. The union of the confi- 
dence masks accepts 73 % of the sky and is denoted U73. It is 
adopted as the default mask for evaluation purposes in this paper. 

In addition to the CMB maps from the full data set, the half- 
ring frequency maps have been processed by each algorithm to 
provide half-ring CMB maps. They are used to provide estimates 
of the instrumental noise contribution to the foreground-cleaned 
maps in the power spectrum analysis (see Sect. 6). The algo- 
rithms were also used to process Monte Carlo simulations: 1000 
realizations of the CMB and 1000 realizations of noise. They 



^ Note that SMICA, being defined in harmonic space, employs a 
smooth filling process inside a small Galactic mask to prevent fore- 
ground residuals from leaking from low to high Galactic latitudes, and 
therefore appears visually diff'erent from the other three solutions in this 
respect; see Appendix D. 



are not used in the analyses presented in this paper, but are used 
by Planck Collaboration XXIII (2013) and Planck Collaboration 
XXIV (2013). 

The beam transfer functions of the foreground-cleaned CMB 
maps have been estimated for each algorithm, as shown in Fig. 4. 
The angular resolution of the NILC, SEVEM, and SMICA maps 
corresponds to a Gaussian beam with FWHM of 5'. The diff'er- 
ence between SEVEM and NILC/SMICA is due to their diff'erent 
treatment of the HEALPix"^ pixel window function (Gorski et al. 
2005). The deviation of NILC beam from a Gaussian shape at 
i > 2800 is caused by the last needlet window (see Appendix B). 
Coimnander-Ruler has a larger beam, because it is defined ex- 
plicitly as a weighted average of frequency maps in pixel space. 
Its resolution is equivalent to a Gaussian beam with FWHM of 
approximately The beam transfer functions have been com- 
puted assuming the best-fit beam transfer function for each fre- 
quency channel, and the uncertainties in the latter have not been 
propagated to these estimates. 

In Fig. 5 we show the standard deviation per pixel among the 
four foreground-cleaned CMB maps downgraded to A^side = 128, 
and in Fig. 6 we show all pairwise diff'erence maps. Typical 
diff'erences at high Galactic latitudes are smaller than 5 jjK. 
Considering the diff'erence maps in more detail, it is clear that 



http : //healpix . sourceforge . net 
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the Commander -Ruler map is the most different from the other 
three, whereas NILC and SMICA are the most similar. This is not 
completely unexpected, because while Commander-Ruler uses 
only frequencies between 30 and 353 GHz in its solution, the 
other three codes additionally include the dust-dominated 545 
and 857 GHz maps. 

This difference in data selection may explain some of 
the coherent structures seen in Fig. 6. In particular, the most 
striking large-scale feature in the difference maps involving 
Commander -Ruler is a large negative band roughly following 
the ecliptic plane. This is where the ZLE (Planck Collaboration 
VI 2013) is brightest. Since the ZLE is also stronger at high fre- 
quencies, having a spectrum close to that of thermal dust, it is 
possible that this pattern may be an imprint of residual ZLE ei- 
ther in the Commander -Ruler map, or in all of the other three 
maps. Both cases are plausible. The Commander-Ruler solu- 
tion may not have enough high-frequency information to distin- 
guish between ZLE and normal thermal dust emission, and, by 
assuming a thermal dust spectrum for the entire high-frequency 
signal at 353 GHz, over- subtracts the ZLE at lower frequencies. 
It is also possible that the other three CMB solutions have posi- 
tive ZLE residuals from extrapolating the high-frequency signal 
model from 857 GHz to the CMB frequencies. Without an ac- 
curate and detailed ZLE model, it is difficult to distinguish be- 
tween these two possibilities. It is of course also possible that 
the true explanation is in fact unrelated to ZLE, and the correla- 
tion with the ecliptic plane is accidental. In either case, it is clear 
that the residuals are small in amplitude, with peak-to-peak val- 
ues typically smaller than lOyuK, of which by far the most is 
contained in a quadrupole aligned with the ecliptic. This pro- 
vides additional evidence that residual ZLE is not important for 
the CMB power spectrum and cosmological parameter estima- 
tion, although some care is warranted when using these maps to 
study the statistics of the very largest angular scales (e.g., Planck 
Collaboration XXIII 2013); checking consistency among all four 
maps for a given application alleviates much of this concern. 

We end this section by showing in Fig. 7 a set of residual 
maps derived by analysing the FFP6 simulation with exactly 
the same analysis approaches as applied to the data. It is evi- 
dent that SMICA produces the map with lowest level of residu- 
als. Considering the morphology in each case, we see that the 
main contaminant for Commander-Ruler is under- subtracted 
free-free emission, while for both NILC and SEVEM it is over- 
subtracted thermal dust emission, and for SMICA it is under- 
subtracted thermal dust emission. However, at high latitudes and 
outside the confidence masks, the residuals are generally below a 
few yuK in amplitude. It is also worth noting that each algorithm 
has been optimized (in terms of model definition, localization 
parameters, etc.) for the data, and the same configuration was 
subsequently used for the FFP6 simulations without further tun- 
ing. The simulations presented here therefore provide a conser- 
vative estimate of the residuals in the data. This is also reflected 
in the fact that the differences between CMB reconstructions for 
the FFP6 simulations are larger than those found in the data. See 
Appendix E for further details. 

6. Power spectrum and cosmological parameters 

In this section we evaluate the foreground-cleaned maps in terms 
of CMB power spectra and cosmological parameters. Our pur- 
pose in doing this is to show that the maps are consistent with 
the high-^ likelihood obtained from the cross-spectrum analy- 
sis of detector set and frequency maps in Planck Collaboration 
XV (2013), and with the cosmological parameters derived from 
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Fig. 8: Angular power spectra of the foreground-cleaned CMB 
maps and half-ring half-difference (HRHD) maps. The spectra 
have been evaluated using the U73 mask apodized with a 30' 
cosine function. 



them in Planck Collaboration XVI (2013). This also establishes 
the consistency between Planck'^ cosmological constraints and 
studies of the large-scale structure and higher order statistics of 
the CMB. 

6.1. Power spectra 

Figure 8 shows the power spectra of the foreground-cleaned 
CMB maps and the corresponding HRHD maps, evaluated using 
the U73 mask with a 30' cosine apodization. The spectra have 
been corrected for the effect of the mask and the beam transfer 
function of each algorithm has been deconvolved. The spectra of 
the HRHD maps give an estimate of the instrumental noise con- 
tribution to the power spectrum of the cleaned map. The correla- 
tions between the HFI half-ring frequency maps are inherited by 
the half-ring CMB maps that use them as input. At small angular 
scales, the CMB solution comes almost entirely from data in the 
HFI channels, and therefore the spectrum of the CMB HRHD 
maps is also biased low. 

At small angular scales, the effective noise levels of NILC, 
SEVEM, and SMICA are very similar, and lower than that of 
Commander -Ruler. The last has larger noise because it operates 
entirely in pixel space and therefore applies the same weights to 
all multipoles. It cannot take advantage of the changing signal- 
to-noise ratio of the frequency channels with angular scale. 

We can estimate the contribution of residual foregrounds to 
the foreground-cleaned CMB maps by making use of the FFP6 
simulations. In addition to processing the simulated frequency 
maps, the maps of the individual input sky components were pro- 
cessed by the algorithms after fixing their parameters or weights 
to the values obtained from the "observed" maps. Figure 9 shows 
the power spectra of the simulated FFP6 components, in this 
case CMB, noise and the sum of the foreground components. 
The top panel shows the spectra computed using the union mask 
derived from the simulation with a 30' cosine apodization. The 
total foreground contribution becomes comparable to the CMB 
signal at ^ 2000. The bottom panel shows the same computed 
with an apodized point source mask applied to the maps (i.e.. 
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Fig. 10: Estimates of the CMB power spectra from the 
foreground-cleaned maps, computed by XF aster. The solid 
lines show the spectra after subtracting the best-fit model of 
residual foregrounds. The vertical dotted line shows the maxi- 
mum multipole = 2000) used in the likelihood for fitting the 
foreground model and cosmological parameters (see Sect. 6.2.2 
for further details). The dashed lines show the spectra before 
residual foreground subtraction. 



Chain Monte Carlo sampler. To speed up this process, we ad- 
ditionally use PICO (Parameters for the Impatient COsmologist, 
Fendt & Wandelt 2008), a tool which interpolates the CMB 
power spectra and matter power spectra as a function of cos- 
mological parameters. 



Fig. 9: Angular power spectra of FFP6 simulated components 
evaluated over the common mask (top) and the common point 
source mask (bottom), both apodized with a 30' cosine func- 
tion. Three components are shown: the CMB (dashed line); noise 
(dot-dashed line); and the sum of all foregrounds (solid line). A 
nonlinear scale is used on the horizontal axis to show all the fea- 
tures of the spectra. 



no diff'use masking, although this mask does removes a large 
part of the Galactic plane). The residual foreground contribution 
is larger at all angular scales, but still it only becomes compa- 
rable to the CMB signal at ^ 1800 in the worst case. For 
both masks, SMICA has the smallest residual foreground con- 
tamination at large angular scales, which is also demonstrated 
in Fig. 7. A more detailed examination of the contribution of the 
individual foreground components to the power spectrum is in 
Appendix E. 

6.2. Likelihood and cosmological parameters 

We estimate the binned power spectra with XFaster (Rocha 
et al. 2009, 2010, 2011) and determine cosmological parame- 
ter constraints using a correlated Gaussian likelihood. Parameter 
constraints are derived using a Metropolis-Hastings Markov 



6.2.1. Model and methods 

We compute the power spectrum for each foreground-cleaned 
map over the multipole range 2 < i < 2500, while parameter 
constraints are derived using only 70 < £ < 2000; as shown 
in Appendix E through simulations, modelling errors become 
non-negligible between f = 2000 and 2500. For parameter es- 
timation, we adopt a standard six-parameter ACDM model, and 
impose an informative Gaussian prior of r = 0.0851 ± 0.014, 
since polarization data are not included in this analysis. 

While the foreground-cleaned maps should have minimal 
contamination from diff'use Galactic emission, they do contain 
significant contamination from unresolved extragalactic sources. 
These contributions are most easily modelled in terms of residual 
power spectra, therefore we marginalize over the corresponding 
parameters at the power spectrum level. To the six ACDM pa- 
rameters, describing the standard cosmology, we add two fore- 
ground parameters, Aps, the amplitude of a Poisson component 
(and hence constant, Q = Aps), and Ad, the amplitude of a clus- 
tered component with shape = l)Q/27r oci^-^ . Both are 
expressed in terms of dXi = 3000 in units of fiYs?. 

The power spectrum calculation is based on the half-ring 
half-sum (HRHS) and HRHD CMB maps (see Sect. 5); the latter 
is used to estimate the noise bias in the power spectra extracted 
from the HRHS maps. From these, we calculate the pseudo- 
spectra, Q and (Hivon et al. 2002), respectively, after apply- 
ing the U73 mask. These are used as inputs to XFaster together 
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Fig. 1 1 : Comparison of cosmological and foreground parameter values estimated from the foreground-cleaned CMB maps for 
^max = 2000 (in red) and those obtained with CamSpec and Plik likelihoods (in blue). The values of the foreground parameters are 
not shown for CamSpec and Plik, since they use a different foreground model. 



To avoid aliasing of power from large to small scales, which 




would add an offset between the signal-plus-noise and noise 
pseudo-spectra at high we use the apodized version of the U73 
mask. The known mismatch in the noise level between the spec- 
tra due to the correlation between the half-ring maps is not ex- 
plicitly corrected. It is left to be absorbed into the two foreground 
parameters. 



Using the pseudo- spectra and XF aster, we then reconstruct 
an estimate of the power spectrum of each foreground-cleaned 
HRHS map, removing the noise bias as estimated from the corre- 
sponding HRHD map. To this end we apply an iterative scheme 
starting from a flat spectrum model. The result is a binned power 
spectrum and the associated Fisher matrix, which are then used 



500 1000 1500 2000 to construct the likelihood, approximated here by a correlated 

Gaussian distribution. 



Fig. 12: Residuals of all map-based best-fit models relative to 

CamSpec best-fit model (assuming a prior on r) for ^max = 2000. - \ p ^ 

^ V to ±^ / max study consistcncy m the low-c range, we nt a two- 

parameter q-n (amplitude-tilt) model relative to the Planck best- 
fit ACDM model on the form, Q = qdlt^^iyoxfCf , using a pixel- 
with the beam transfer functions provided by each method (see space likelihood for maps smoothed to 6° FWHM; see Planck 
Fig. 4). Collaboration XV (2013) for further algorithmic details. 
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Fig. 13: Low-^ power spectrum amplitude and tilt constraints 
measured relative to the best-fit Planck ACDM model derived 
from foreground-cleaned CMB maps smoothed to 6° FWHM. 
The cross shows the best-fit model {q,n) = (1,0). 



6.2.2. Results 

We perform the power spectrum and parameter estimation anal- 
ysis for both the data and the FFP6 simulations described in 
Sect. 4. The results for the latter are given in Appendix E. 

Figure 10 shows estimates of the angular power spectrum 
for each foreground-cleaned map, with the uncertainties given 
by the Fisher matrix. The parameter summary given in Fig. 1 1 
shows the parameter constraints derived using multipoles be- 
tween £ = 70 and 2000, and compares these to results obtained 
with the CamSpec and Plik likelihoods (Planck Collaboration 
XV 2013). 

Diff'erences in the power spectra at high £ are mostly ab- 
sorbed by the two-parameter foreground model, rendering con- 
sistent cosmological parameters. For example, the high-^ power 
excess seen in the Coininander- Ruler map is well-fitted in terms 
of residual point sources, which makes intuitive sense, consid- 
ering the lower angular resolution of this map (see Sect. 5). 
However, the ACDM parameter uncertainties derived from the 
four codes are very consistent. This indicates that most of the 
cosmological information content above £ > 1500 is degenerate 
with the extragalactic foreground model, and a more sophisti- 
cated foreground treatment is required in order to recover sig- 
nificant cosmological parameter constraints from these scales. 
Beyond this, deviations among cosmological parameters are 
small and within 1 cr for all methods and most of the parame- 
ters. Further, the parameters derived from the four foreground- 
cleaned CMB maps are in good agreement with those obtained 
by CamSpec and Plik using cross-spectra; departures are well 
within 1 cr for most parameters. 

Inspecting the diff'erences between the best-fit models de- 
rived from the four foreground-cleaned maps and from CamSpec 
plotted in Fig. 12, we find that the relative residuals are within 
40 yuK^ for all multipole ranges, and smaller than 20 yuK^ at high 
£. This can be compared to the corresponding residuals for the 
FFP6 simulation shown in Appendix E. 

The likelihood used for this analysis does not take into ac- 
count some systematic eff'ects that will aff'ect our foreground- 
cleaned CMB maps, such as relative calibration uncertainties 



between the frequency channel maps used to construct them, 
or their beam uncertainties. These eff'ects are accounted for in 
the likelihoods in Planck Collaboration XV (2013). We have 
also adopted a very simple two-parameter model for the resid- 
ual extragalactic foregrounds. Despite these limitations, the four 
CMB maps yield cosmological parameters in agreement with 
the cross-spectrum based likelihoods for a basic six-parameter 
ACDM model. Thus we can be confident that the CMB maps 
are consistent with the power spectrum analysis. 

Before concluding this section, we show in Fig. 13 the results 
from a two-parameter fit of an amplitude-tilt model to each of 
the four foreground-cleaned maps, downgraded to 6° and repix- 
elized at an A^side = 32 grid. Clearly, the maps are virtually identi- 
cal on large angular scales measured relative to cosmic variance, 
with any diff'erences being smaller than 0.1 cr in terms of cosmo- 
logical parameters. However, it is worth noting that the best-fit 
model, (q,n) = (1,0), is in some tension with the low-^ spec- 
trum, at about 1.7 cr in this plot. The same tension between large 
and small angular scales is observed in Planck Collaboration XV 
(2013) and Planck Collaboration XVI (2013) with higher statis- 
tical significance using the full Planck likelihood. Irrespective of 
physical interpretation, the calculations presented here demon- 
strate that these low-^ features are robust with respect to compo- 
nent separation techniques. 



7. Higher-order statistics 

The foreground-cleaned CMB maps presented in this pa- 
per are used as inputs for most Planck analyses of higher- 
order statistics, including non-Gaussianity studies (Planck 
Collaboration XXIV 2013), studies of statistical isotropy 
(Planck Collaboration XXIII 2013), gravitational lensing by 
large-scale structure (Planck Collaboration XVII 2013), and of 
the integrated Sachs- Wolfe effect (Planck Collaboration XXIV 
2013). In this section we provide a summary of the non- 
Gaussianity and gravitational lensing results. 

7.1. Non-Gaussianity 

Primordial non-Gaussianity is typically constrained in terms of 
the amplitude, f^l^^, of the quadratic corrections to the gravi- 
tational potential, as well as by means of the three-point cor- 
relation function based on diff'erent triangle configurations. The 
results from these calculations for the foreground-cleaned CMB 
maps are presented in Planck Collaboration XXIV (2013). After 
subtraction of the lensing-ISW correlation contribution, the final 
result is f^^^^ = 2.7+5.8, as estimated from the SMICA map using 
the KSW bispectrum estimator (Komatsu et al. 2005), consis- 
tent within 1 cr with results from other methods and foreground- 
cleaned maps. 

Uncertainties are evaluated by means of the FFP6 simula- 
tions, and potential biases are studied using both Gaussian and 
non-Gaussian CMB realizations. In particular, when a detectable 
level of primordial non-Gaussianity (f^l^^ = 20.4075) is in- 
jected into the FFP6 simulations, each foreground-cleaned map 
yielded a positive detection within 2 cr of the expected value, 
recovering values of f^l^^ = 8.8 ± 8.6, 19.0 ± 7.5, 11.1+7.6 
and 19.7 + 7.4 for Commander-Ruler, NILC, SEVEM, SMICA, re- 
spectively. We see that NILC and SMICA demonstrate the best 
recovery of the injected non-Gaussianity, and we favoured the 
latter for non-Gaussian studies for its faster performance over 
NILC. The foreground-cleaned CMB maps presented in this pa- 
per do not provide significant evidence of a non-zero value of 
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Fig. 14: Lensing power spectrum estimates from FFP6 simula- 
tions using an apodized mask covering /sky 2 - 0.70 of the sky. 



f^^^^, and realistic simulations show that the component sep- 
aration methods do not suppress real non-Gaussian signatures 
within expected uncertainties. The implications of these results 
in terms of early Universe physics are discussed in the relevant 
papers (Planck Collaboration XXIV 2013; Planck Collaboration 
XXII 2013). 



7.2. Gravitational lensing by large-scale structure 

Gravitational lensing by the intervening matter imprints a non- 
Gaussian signature in the CMB, which allows the reconstruction 
of the gravitational potential integrated along the line of sight to 
the last scattering surface. In Planck Collaboration XVII (2013), 
this effect has been detected at a high significance level (greater 
than 25 cr) using the Planck temperature maps. Specifically, the 
lensing induced correlations between the total intensity and its 
gradients have been used to reconstruct a nearly full sky map of 
the lensing potential 0, which has been used for further studies 
on Planck data, including the detection of a non-zero correla- 
tion with the ISW (Planck Collaboration XXIV 2013; Planck 
Collaboration XIX 2013) and other tracers of large-scale struc- 
ture (notably, significant correlation with the CIB is reported in 
Planck Collaboration XVIII 2013), as well as the estimate of 
the power spectrum of the lensing potential and the associated 
likelihood. The latter was constructed using a simple minimum 
variance combination of the 143 and 217 GHz maps on about 
70 % of the sky, as well as subtracting dust contamination using 
the 857 GHz Planck channel as a template (Planck Collaboration 
XVII 2013). These lensing results have improved the cosmolog- 
ical constraints from Planck (Planck Collaboration XVI 2013). 

The foreground-cleaned CMB maps described in Sect. 5 
were used to perform a lensing extraction on a larger sky frac- 



tion, reaching about 87 % of the sky. We found the lensing power 
spectrum to be in good agreement with the one obtained using 
the minimum variance combination, i.e., the signal agrees within 
1 cr in the majority of the angular domain bins, and is character- 
ized by an equivalent uncertainty. The foreground-cleaned maps 
were further exploited on the baseline 70 % sky fraction for as- 
sessing the robustness of the main reconstruction against the 
foreground contamination (Planck Collaboration XVII 2013). 

We show that the component separation algorithms pre- 
sented in this paper do not bias the lensing reconstruction in 
the case of the large sky fraction considered here. We consider 
FFP6 simulations including noise and lensed CMB signal, prop- 
agated through each of the component separation algorithms de- 
scribed in Sect. 3. We perform a lensing potential reconstruction 
in the pixel domain based on the CMB maps processed by the 
four component separation methods using the metis algorithm 
described in Planck Collaboration XVII (2013). This method 
uses the quadratic estimator presented in Okamoto & Hu (2003), 
which corrects for the mean-field bias caused by extra sources of 
statistical anisotropy in addition to the CMB. 

For each method, we combine the masks of CO regions, 
nearby galaxies and compact objects as defined in Planck 
Collaboration XVII (2013), with the CG90 mask described in 
Sect. 4. This procedure results in masks with sky fractions 
/sky = 0.836,0.851,0.850,0.846 for Coimnander -Ruler, NILC 
SEVEM, and SMICA, respectively. 

We estimate the lensing potential power spectrum, C^^, 
following the methodology described in Planck Collaboration 
XVII (2013). It consists of a pseudo-Q estimate based on 
a highly-apodized version of the lensing potential reconstruc- 
tion, which has an eff'ective available sky fraction /sky,2 = 
0.648,0.690,0.686,0.683 for Coimnander -Ruler, NILC, SEVEM 
and SMICA, respectively. The band-power reconstructions in 17 
bins in the range 2 < i < 1025 are plotted in Fig. 14, as well as 
the residuals relative to the theoretical lens power spectrum. All 
algorithms achieved an unbiased estimation of the underlying 
lensing power spectrum, with = 10.58, 17.34, 18.54, 15.30, 
for Coimnander -Ruler, NILC, SEVEM, and SMICA respectively, 
with 17 degrees of freedom. The associated probability-to- 
exceed (PTE) values are 83 %, 36 %, 29 %, 50 %. 

The power spectrum estimates are in remarkable agreement 
with each other. However, the Commander -Ruler solution has 
significantly larger uncertainties, as expected from its lower 
signal-to-noise ratio to lensing due to its larger beam. These 
results on simulated foreground-cleaned CMB maps demon- 
strate that the component separation algorithms do not alter the 
lensing signal, and this provides a strategy for achieving a ro- 
bust lensing reconstruction on the largest possible sky cover- 
age. The foreground-cleaned maps have been used in Planck 
Collaboration XVII (2013) to obtain lensing potential estimates 
on 87% of the sky. 



8. Foreground components 

In this section we consider the diff'use Galactic components, and 
present full-sky maps of thermal dust and CO emission, as well 
as a single low-frequency component map representing the sum 
of synchrotron, AME, and free-free emission. Our all-sky CO 
map is a "type 3" product as presented in Planck Collaboration 
XIII (2013). To assess the accuracy of these maps, we once again 
take advantage of the FFP6 simulation. The Commander -Ruler 
method used in the following is described in Appendix A and 
consists of a standard parametric Bayesian MCMC analysis at 
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Fig. 15: Posterior mean foreground amplitude maps derived 
from the low-resolution analysis. From top to bottom are shown 
the low-frequency, CO and thermal dust emission maps. 



low angular resolution, followed by a generalized least-squares 
solution for component amplitudes at high resolution. 

8.1. Data selection and processing 

We only use the seven lowest Planck frequencies, from 30 to 
353 GHz. The two highest channels have significantly differ- 
ent systematic properties than the lower frequency bands, for 
instance concerning calibration, ZLE, and noise correlations. 
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Fig. 16: Posterior mean spectral parameter maps derived from 
the low-resolution analysis. The top panel shows the power law 
index of the low-frequency component, and the bottom panel 
shows the emissivity index of the one-component thermal dust 
model. Note that the systematic error due to monopole and 
dipole uncertainties is significant for the dust emissivity in re- 
gions with a low thermal dust amplitude. 



and they are more relevant to thermal dust and CIB studies 
than to the present CMB analysis. Studies of specific fore- 
grounds (CO, thermal dust, CIB etc.) using all Planck frequen- 
cies as well as ancillary data are discussed in companion pa- 
pers Planck Collaboration VI (2013), Planck Collaboration XXI 
(2013), Planck Collaboration (2013b), and additional future pub- 
lications will consider extensions to AME, synchrotron and free- 
free emission. 

In order to obtain unbiased estimates of the spectral param- 
eters across all frequency bands, each map is downgraded from 
its native resolution to a common angular resolution of 40' and 
repixelized at A^side = 256, a limit imposed by the LFI 30 GHz 
channel. Once the spectral indices have been determined, we re- 
estimate the component amplitudes at native Planck resolution 
(see Appendix A). 

Although the smoothing operation introduces noise correla- 
tions between pixels, we model the noise of the smoothed maps 
as uncorrected white noise with an eff'ective standard deviation, 
cr(p), for each pixel p. This approximation does not bias the fi- 
nal solution, because the analysis is performed independently 
for each pixel. However, it is important to note that correlations 
between pixels are not taken into account in this analysis. The ef- 



13 



Planck Collaboration: Planck 2013 results. XIL Component separation 




14 



Fig. 17: P^r pixel for the joint CMB and foreground analysis. 
The expected value for an acceptable fit is 7, corresponding to 
the number of frequency bands used in this analysis. The pixels 
with high values can be classified into two types, due to either 
modelling errors (i.e., high residuals in the Galactic plane) or to 
un-modelled correlated noise (i.e., stripes crossing through low 
dust emission regions). 



fective noise uncertainty, cr{p), is estimated using realistic noise 
simulations downgraded in the same way as the data. The mea- 
sured instrumental bandpasses are taken into account by inte- 
grating the emission laws over the bandpass for each component 
at each Monte Carlo step in the analysis. 

The monopole (zero-point) of each frequency map is not 
constrained by Planck, but is rather determined by post- 
processing, and associated with a non-negligible uncertainty (see 
Table 5 of Planck Collaboration I 2013). In addition, each fre- 
quency map includes a significant monopole contribution from 
isotropic extragalactic sources and CIB fluctuations not traced 
by local Galactic structure, ranging from less than about 10- 
20 yuK at 70 GHz to several hundreds of yuK at 353 GHz. Finally, 
the efl'ective dipole in each map is associated with significant 
uncertainty due to the large kinematic CMB dipole. In order to 
prevent these eff'ects from introducing modelling errors during 
component separation, they must be fit either prior to or jointly 
with the Galactic parameters. Unfortunately, when allowing free 
spectral parameters per pixel, there is a near-perfect degeneracy 
among the off'sets, the foreground amplitudes and the spectral 
indices, and in order to break this degeneracy, it is necessary to 
reduce the number of spectral degrees of freedom. 

We adopt the method described by Wehus et al. (2013) 
for this purpose, which has the additional advantage of mak- 
ing minimal assumptions about the foreground spectra. In short, 
this method uses linear regression between data from CMB- 
subtracted maps evaluated on pixels falling within each large 
^side = 8 pixel to estimate the relative off'sets, mi and m2, be- 
tween any two maps at each position on the sky. Each regression 
provides a constraint of the form m\ = am2 + b, where a and b 
are the slope and off'set, respectively, and where each value of m/ 
consists of the sum of both a monopole and a dipole term evalu- 
ated at that position. The individual monopoles and dipoles can 
then be reconstructed by measuring a and b in diff'erent regions 
of the sky, exploiting spatial variations in spectral indices, and 
solving jointly for two monopoles and dipoles, including con- 
straints from all positions. To minimize degeneracies, a positiv- 
ity prior is imposed on the fit, such that statistically significant 
negative pixels are heavily penalized. For 44 and 70 GHz, we 



retain the dipole values determined during the mapmaking pro- 
cess, and do not attempt to fit them. 

The resulting complete set of monopole and dipole values is 
listed in Table 2. As a cross-check, we performed a dedicated 
Coimnander run in which we fitted for the dipole at 353 GHz, to- 
gether with the foreground amplitudes and spectral indices, and 
only found sub-yuK diff'erences. This channel is by far the most 
problematic in our data set in terms of off'set determination, be- 
cause of the very bright dust emission at this frequency. As a 
result, there is a large relative uncertainty between the zero-level 
of the dust amplitude map and the 353 GHz channel off'set not 
accounted for in the following analyses. However, the sum of 
the two terms is well determined, and a potential error in ei- 
ther therefore does not compromise the quality of the other sig- 
nal components (e.g., CMB and low-frequency components). A 
similar comment applies between the off'set at 30 GHz and the 
zero-level of the low-frequency component, although at a signif- 
icantly lower level. 

8.2. Component models and priors 

Our model for the low-resolution CMB analysis includes four in- 
dependent physical components: CMB; "low-frequency" emis- 
sion; CO emission; and thermal dust emission. It can be written 
schematically in the form 



Sy(p) = Acmb(p) + Aif(p) 

\^o,if/ 

- 1 / y \^'^^^^+^ 

+ Aco(p)fv,co+A^(p)— — , (1) 

^kT^ip) _ 1 \^0,d/ 

where Ai(p) denotes the signal amplitude for component / at 
pixel p, vqj is the reference frequency for each component, and 
y refers to frequency. (Note that for readability, integration over 
bandpass, as well as unit conversions between antenna, flux den- 
sity and thermodynamic units, is suppressed in this expression.) 
Thus, each component is modelled with a simple frequency 
spectrum parameterized in terms of an amplitude and a small 
set of free spectral parameters (a power-law index for the low- 
frequency component, and an emissivity index and temperature 
for the thermal dust component); no spatial priors are imposed. 
One goal of the present analysis is to understand how well this 
simple model captures the sky signal in terms of efl'ective com- 
ponents over the considered frequency range, and we exploit the 
FFP6 simulation (see Sect. 4) for this purpose. 

In order to take into account the efl'ect of bandpass integra- 
tion, each term in the above model is evaluated as an integral 
over the bandpass as described in Sect. 3 of Planck Collaboration 
IX (2013), and converted internally to thermodynamic units. 
Accordingly, the reference frequencies in Eq. 1 are computed 
as efl'ective integrals over the bandpass, such that the amplitude 
map, Ai(p), corresponds to the foreground map observed by the 
reference detector, i.e., after taking into account the bandpass. 
In order to minimize degeneracies between the difl'erent signal 
components, the reference band for a given component is set to 
the frequency at which its relative signal-to-noise ratio is maxi- 
mized. 

The foreground model defined in Eq. 1 is motivated by prior 
knowledge about the foreground composition over the CMB fre- 
quencies as outlined in Sect. 2, as is our choice of priors. In 
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Fig. 18: Comparison of the high-resolution Ruler (top) and low-resolution Coininander (bottom) amplitude maps for a particularly 
strong CO complex near the Fan region; the maps are centred on Galactic coordinates (/,/?) = (110°, 15°), and the grid spacing is 
5°. Columns show, left to right: the CMB amplitude; the CO amplitude at 100 GHz and the thermal dust amplitude at 353 GHz. 

Table 2: Estimated monopoles and dipoles in Galactic coordinates, all measured in thermodynamic yuK. Errors are estimated by 
bootstrapping, and do not account for correlated errors across frequencies. In particular, the 353 GHz monopole uncertainty is 
dominated by systematic errors not included in these estimates. Note that the dipoles at 44 and 70 GHz are fixed at the values 
determined in the mapmaking. 



Frequency Monopole X dipole Y dipole Z dipole 

[GHz] \m im im 

30 8 + 2 -4 + 3 -6 + 2 6 + 1 

44 2 + 1 + + + 

70 15 + 1 + + + 

100 15 + 1 2 + 1 5 + 1 -5 + 1 

143 33 + 1 2 + 1 7 + 1 -6 + 1 

217 86 + 1 2 + 1 11+2 -10 + 2 

353 414 + 4 11 + 10 52 + 12 -37 + 8 



addition to the Jeff'reys prior^ (Eriksen et al. 2008), we adopt 
Gaussian priors on all spectral parameters with centre values and 
widths attempting to strike a balance between prior knowledge 
and allowing the data to find the optimal solution. Where needed, 
we have also run dedicated analyses, either including particular 
high signal-to-noise ratio subsets of the data or using a lower res- 

^ The purpose of the Jeffreys prior is to normalize the parameter vol- 
ume relative to the likelihood, such that the likelihood becomes so- 
called "data- translated", i.e., invariant under re-parameterizations. 



olution parameterization to increase the eff'ective signal-to-noise 
in order to inform our prior choices. We now consider each fore- 
ground component in turn, and note in passing that the CMB 
component, by virtue of being a blackbody signal, is given by a 
constant in thermodynamic temperature units. 

We approximate the low-frequency component by a straight 
power law in antenna temperature with a free spectral index per 
pixel, and adopt a prior of p = -3 + 0.3 (this is the index in 
terms of brightness temperature). This choice is determined by 
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noting that the prior is in practice only relevant at high Galactic 
latitudes where the signal-to-noise ratio is low and the domi- 
nant foreground component is expected to be synchrotron emis- 
sion; in the signal-dominated and low-latitude AME and free- 
free regions, the data are sufficiently strong to render the prior 
irrelevant. For validation purposes, we have also considered mi- 
nor variations around this prior, such sls JS = -2.9 ± 0.3 and 
y5 = -3.05 ±0.2, finding only small diff'erences in the final so- 
lutions. The reference band for the low-frequency component is 
set to 30 GHz, where the low-frequency foreground signal peaks. 
The final low-frequency amplitude map is provided in units of 
thermodynamic microkelvin. 

The CO emission is modelled in terms of a single line ratio 
for each frequency. Specifically, the CO amplitude is normalized 
to the 100 GHz band, and defined in units of yuKkm s"^ (Planck 
Collaboration XIII 2013). The amplitude at other frequencies is 
determined by a single multiplicative factor relative to this, with 
a numerical value of 0.595 at 217 GHz and 0.297 at 353 GHz; 
all other frequencies are set to zero. These values are obtained 
from a dedicated CO analysis that includes only high signal-to- 
noise ratio CO regions covering a total of 0.5% of the sky. The 
derived values are in good agreement with those presented by 
Planck Collaboration XIII (2013). 

Thermal dust emission is modelled by a one-component 
modified blackbody emission law with a free emissivity spectral 
index, y^d, and dust temperature, Td, per pixel. However, since 
we only include frequencies below 353 GHz, the dust temper- 
ature is largely unconstrained in our fits, and we therefore im- 
pose a tight prior around the commonly accepted mean value of 
Td = 18 ± 0.05 K. The only reason we do not fix it completely 
to 18 K is to allow for modelling errors near the Galactic centre. 
The dust emissivity prior is set toy^d = 1-5 ± 0.3, where the mean 
is determined by a dedicated run fitting for a single best-fit value 
for the high-latitude sky, where the prior is relevant. The refer- 
ence band for the thermal dust component is 353 GHz, and the 
final map is provided in units of megajansky per steradian. 

8.3. Results and validation 

The output of the Bayesian component separation algorithm 
is a set of samples drawn from the joint posterior distribution 
of the model parameters, as opposed to a single well-defined 
value for each. For convenience, we summarize this distribution 
in terms of posterior mean and standard deviation maps, com- 
puted over the sample set, after rejecting a short burn-in phase. 
The goodness-of-fit is monitored in terms of the x'^ P^r pixel. 
Although convenient, it is, however, important to note that this 
description does not provide a comprehensive statistical repre- 
sentation of the full posterior distribution, which is intrinsically 
non-Gaussian. One should be careful about making inferences in 
the low signal-to-noise regime based on this simplified descrip- 
tion. 

The low-resolution Coimnander posterior mean amplitude 
maps are shown in Fig. 15 for the low-frequency, CO, and ther- 
mal dust components, and the spectral index maps in Fig. 16. 
The associated map is plotted in Fig. 17. Note that because 
we are sampling from the posterior instead of searching for the 
maximum-likelihood point, the expected number of degrees of 
freedom is equal to A^band = 7 in this plot, not A^band - Np^r- 

Figure 18 compares the high-resolution Ruler solution to 
the low-resolution Coimnander solution for CMB, CO and ther- 
mal dust on a particularly strong CO complex near the Fan re- 
gion, centred on Galactic coordinates (/,/?) = (110°, 15°). 




-30 AiK 30 

(a) Low-frequency component residual at 30 GHz 
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(b) CO residual at 100 GHz 
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(c) Thermal dust residual at 353 GHz 



Fig. 19: Amplitude residual maps, Aout - ^in, computed blindly 
from the FFP6 simulation. The panels show (from top to bot- 
tom) the low-frequency residual at 30 GHz, the CO residual at 
100 GHz and the thermal dust residual at 353 GHz. All units are 
thermodynamic yuK. The white lines indicate the boundary of the 
Coimnander likelihood analysis mask, removing 13% of the sky. 



Several features can be seen here, foremost of which is that 
the Galactic plane is strikingly obvious, with;^^ values exceed- 
ing 10^ for seven degrees-of-freedom in a few pixels. This is 
not surprising, given the very simplified model at low frequen- 
cies (i.e., a single power law accounting for AME, synchrotron, 
and free-free emission), as well as the assumption of a nearly 
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constant dust temperature of 18 K. Second, there is an extended 
region of moderately high;^^ roughly aligned with a great circle 
going through the Ecliptic South Pole, indicating the presence 
of correlated noise in the scanning rings not accounted for in our 
white noise model. 

Based on these, and other considerations, it is clear that parts 
of the sky must be masked before proceeding to CMB power 
spectrum and likelihood analyses. This masking process is dis- 
cussed at greater length in Planck Collaboration XV (2013), and 
results in different masks for specific applications. The goal of 
our present discussion is to evaluate the adequacy of the mask 
adopted for low-^ likelihood analysis (L87), which is based on 
the fits presented here. This mask removes 13% of the sky, and 
is derived from a combination of and component amplitude 
thresholding. 

For validation purposes, we analyse the simulations de- 
scribed in Section 4 in the same way as the real data, includ- 
ing monopole and dipole determination, CO line ratio estima- 
tion and spectral index estimation. Individual component maps 
at each observed frequency are available from the simulation 
process, and used for direct comparison with the reconstructed 
products. 

In Fig. 19 we show the diff'erences between the recovered 
and input component maps at their respective reference frequen- 
cies. The boundary of the 13% Coimnander mask is traced by the 
white contours, and a best-fit monopole and dipole have been 
subtracted from each diff'erence map. All diff'erence maps are 
shown in units of thermodynamic yuK. The top panel of Fig. 20 
gives the error histograms outside the masked region for each 
component, normalized to the respective estimated standard de- 
viation; if the recovered solution has both correct mean and stan- 
dard deviation, these histograms should match a Gaussian distri- 
bution with zero mean and unit variance, indicated by the dashed 
black line. Conversely, a significant bias would be visible as a 
horizontal shift in this plot, while under-estimation of the errors 
would result in too wide a distribution and vice versa. The bot- 
tom panel shows the fractional error (i.e., the error divided by 
the true input value) for all pixels with signal above 5cr\ the 
fractional error is not a useful quantity for noisier signals. 

The diff'erence maps in Fig. 19 display significant errors in 
the Galactic plane. For the low-frequency component, the resid- 
uals are dominated by free-free emission, while for thermal dust 
the dominant contaminant is CO emission. However, outside the 
mask the residuals are small, and, at least for the low-frequency 
and CO components, the spatial characteristics appear similar to 
instrumental noise. This is more clear in the histograms shown 
in the top panel of Fig. 20; the mean and standard deviations are 
^if = 0.01 + 1.12, ^co = 0.00+0.87, and^td = 0.00+2.01, respec- 
tively, for the low-frequency, CO and thermal dust components. 
There is no evidence of bias outside the mask in any component, 
and the error estimates are accurate to 12 and 13 % for the low- 
frequency and CO components. Note, though, that the estimated 
error for the CO component is actually larger than the true un- 
certainty, suggesting that the white noise approximation for the 
100 GHz channel overestimates the true noise. This can occur if 
the correlated instrumental noise is important in regions where 
there is no significant CO emission. Locally re-scaling the white 
noise to account for spatially varying correlated noise would cor- 
rect this effect. 

For the thermal dust component, on the other hand, the er- 
ror is underestimated by a factor of 2. The explanation for this 
is most easily seen from the lower panel of Fig. 19. This map 
is dominated by isotropic CIB fluctuations, rather than instru- 
mental noise. Because these fluctuations have a slightly differ- 
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Fig. 20: Error validation for component amplitudes, evaluated 
from the FFP6 simulation. The upper panel shows histograms 
of the normalized errors 6 = (Aout - Ain)/<^out for the three 
foreground components and including all pixels outside the 
Coimnander likelihood analysis mask. The lower panel shows 
histograms of the fractional error / = (Aout - Ain)/Ain for pix- 
els with a foreground detection level above 5 cr. No evidence of 
significant bias is observed for any component, and the uncer- 
tainty estimates for the low-frequency and CO components are 
accurate to about 12%; the thermal dust uncertainty is under- 
estimated by a factor of 2 due to the presence of unmodelled 
fluctuations. 



ent spectrum than the dominant Galactic dust emission, and the 
model does not account for a separate CIB component, the error 
on the Galactic component is underestimated. When using the 
Galactic map presented here for detailed analysis near the noise 
limit, taking into account these residual fluctuations is essential, 
and the effective noise per pixel should be increased by a factor 
of 2. 

As clearly seen in Fig. 19, the residuals inside the mask are 
highly significant in a strict statistical sense. However, as seen in 
the bottom panel of Fig. 20, they are relatively small in terms of 
fractional errors. Specifically, the three histograms have means 
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Fig. 21: Validation of spectral parameters for low-frequency 
foregrounds, thermal dust, and CO emission, evaluated from the 
FFP6 simulation. Each histogram shows the error distribution 
at the two leading sub-dominant frequencies in the form of the 
normalized errors 6 = (Aout(v) - Ain(y))/(Tout(v) for all pixels 
outside the Commander likelihood analysis mask, where Aout(v) 
is the predicted foreground amplitude at frequency v given the 
estimated amplitude and spectral parameters, and crout(v) is the 
corresponding standard deviation computed over the sample set. 



and standard deviations of fy = 0.00 ± 0.10, fco = -0.03 ± 0.10, 
and /td = 0.00 ± 0.06, respectively, for the low-frequency, CO 
and thermal dust components. The largest bias is observed for 
the CO component, for which the absolute amplitude is biased 
by 3 %. The bias in the low-frequency and thermal dust com- 
ponents is negligible, and the fractional uncertainties are 10 and 
6 %, respectively. This confirms that approximating the sum of 
the three low-frequency components by a single power-law over 
the Planck frequency bands is reasonable; if modelling errors 
dominated, one would expect to see a significant bias in the re- 
sulting amplitude. 

In order to validate the spectral parameters, we show in 
Fig. 21 histograms of the normalized residuals for each fore- 
ground component evaluated at its two leading sub-dominant 
frequencies (i.e., at 44 and 70 GHz for the low-frequency com- 
ponent; at 217 and 353 GHz for the CO component; and at 143 
and 217 GHz for the thermal dust component). The means and 
standard deviations of these distributions are: ^if(44GHz) = 
-0.41 ± 1.98 and ^if(70GHz) = -0.34 ± 2.04 for the low- 
frequency component; ^co(217GHz) = 0.10 ± 0.84 and 
Jco(353GHz) = 0.51 ± 1.00 for the CO component; and 
(5td(143GHz) = -0.02 ± 1.53 and ^td(217GHz) = -0.13 ± 1.87 
for the thermal dust component. As expected, the eff'ect of mod- 
elling errors is more significant at the sub-dominant frequen- 
cies than at the pivot frequencies, when measured in terms of 
statistical uncertainties, since the foreground signal is weaker 
and the confusion with the other components relatively larger. 
Nevertheless, we see that the absolute bias is at most 0.5 cr for 
the CO component at 353 GHz, while the thermal dust bias is 
negligible even at 143 GHz. The estimated uncertainties are gen- 
erally underestimated by up to a factor of two due to these mod- 
elling errors. 
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Fig. 22: Top: Diff'erence map between the estimated thermal 
dust amplitude at 353 GHz derived by Planck Collaboration 
(2013b) and the low-resolution dust map presented here, both 
smoothed to 40'. A monopole and dipole term outside the 
Commander mask has been removed. The former includes only 
high-frequency observations (the Planck 353, 545 and 857 GHz 
channels and observations at 100 yum), while the one derived 
here only uses low-frequency data (Planck 30-353 GHz). Note 
that the colour scale ranges between -0.02 and 0.02 MJy sr"^ in 
this plot, whereas it ranges from to 2.5 MJy sr"^ in the bottom 
panel of Fig. 15. Masked regions indicate pixels with a CO am- 
plitude at 100 GHz larger than lyuKkms"^. Bottom: T-T plot 
between the same two maps. 



Finally, the efficiency of the adopted foreground model for 
CMB analysis is quantified in Appendix E in terms of power 
spectrum residuals and cosmological parameter estimation. 

To summarize, we find that the simplified model, defined by 
Eq. 1, provides a good fit to the realistic FFP6 simulation for 
most of the sky. Absolute residuals are small, and the amplitude 
uncertainty estimates are accurate to around 12 %, except for the 
thermal dust component for which unmodelled CIB fluctuations 
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are important. Further, we find that the real Planck data behave 
both qualitatively and quantitatively very similarly to the FFP6 
simulation, suggesting that this approach also performs well on 
the real sky. 

8.4. Interpretation and comparison with other results 

The maps shown in Figs. 15 and 16 provide a succinct summary 
of the average foreground properties over the Planck frequency 
range. We now consider their physical interpretation and com- 
pare them to products from alternative methods. 

First, the top panel of Fig. 22 shows a diff'erence map be- 
tween the dust map at 353 GHz derived in the present paper 
and one determined from only the three highest Planck frequen- 
cies and the 100 jim IRIS map by Planck Collaboration (2013b), 
shown in flux density units, after removing a small monopole 
and dipole difl'erence. This map is to be compared to the corre- 
sponding dust amplitude map in the bottom panel of Figure 15; 
note that the colour range varies from -0.02 to 0.02 MJy sr"^ in 
the difl'erence plot, and between and 2.5 MJy sr"^ in the ampli- 
tude plot. Thus, despite the very difl'erent data sets and methods, 
we see that the two reconstructions agree to very high accuracy 
outside the Galactic plane. Inside the Galactic plane, the difl'er- 
ences are dominated by residuals due to difl'erent CO modelling, 
seen as solid blue colours in Fig. 22; however, even in this re- 
gion the difl'erences are smaller than 5 % of the amplitude. The 
bottom panel shows a corresponding T-T plot between the two 
maps, excluding any pixel for which the Coimnander CO ampli- 
tude at 100 GHz is larger than 1 yuKkm s"^ The two maps agree 
to 0.2% in terms of best-fit amplitude. 

From Fig. 16, we see that the dust emissivity ranges between 
1.3 and 1.7 for most of the sky; considering only the pixels with 
a posterior distribution width that is a third of the prior width 
(i.e., (T(ySd) < 0-1), we find a mean value of 1.49. The two ex- 
ceptions are a large region of shallow indices northeast of the 
Galactic centre, and steep indices near the Galactic plane. The 
former region corresponds to a part of the sky with low dust 
emission, where we expect the spectral index to be sensitive to 
both monopole and dipole residuals, as well as instrumental sys- 
tematics, such as correlated 1 // noise. The latter appears to be 
particularly pertinent here because the shallow index region at 
least partially traces the Planck scanning strategy; as a result, the 
systematic error on the spectral index in this region is consider- 
able. The main systematic uncertainty connected to the region 
of steep indices around the Galactic plane is confusion with CO 
emission. 

The CO map shown in Figs. 15 and 18 is discussed in greater 
detail in Planck Collaboration XIII (2013). A distinct advan- 
tage of this particular solution over available alternatives is its 
high signal-to-noise ratio per pixel, which is achieved by reduc- 
ing all information into a single value per pixel. Consequently, 
this map serves as a unique tool for follow-up CO observations. 
However, the assumption of a constant line ratio over the full sky 
may lead to a significant systematic uncertainty on CO ampli- 
tude per pixel. This possibility is investigated in a forthcoming 
work (Planck Collaboration 2013a). 

Finally, the spectral index map for the low-frequency com- 
ponent shown in Fig. 16 can be used to determine the domi- 
nant low-frequency component (synchrotron, free-free or AME) 
as a function of position on the sky. To illustrate this connec- 
tion, we once again take advantage of the FFP6 simulation for 
which we know the amplitude of each low-frequency compo- 
nent per pixel. In the top panel of Fig. 23, we use this informa- 
tion to make a "dominant component map"; dark blue indicates 
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Fig. 23: Top: Dominant foreground component per pixel at 
30 GHz in the FFP6 simulation. Dark blue indicates that syn- 
chrotron emission is the strongest component at 30 GHz, light 
blue indicates that free-free dominates, and orange indicates 
that spinning dust (AME) is the strongest component. Bottom: 
The recovered low-frequency power-law index derived from the 
same simulation. 



that synchrotron emission is strongest at a given pixel, light blue 
that free-free is strongest, and orange that spinning dust (AME) 
dominates. In the bottom panel, we show our derived power-law 
index map from the same simulation. As expected, the corre- 
spondence between the power-law index and the dominant com- 
ponent is very strong, implying that the spectral index map can 
be used to trace the individual components. In particular, we 
see that an index below about -3.3 reflects the presence of a 
component consistent with a spinning dust model peaking below 
20 GHz over the Planck frequency range, while an index higher 
than around -2.3 signals the importance of free-free emission. 
Intermediate values typically indicate synchrotron emission, al- 
though it should be noted that the signal-to-noise ratio at very 
high latitudes is low and the results are therefore prior-driven in 
these regions. 

Returning to the spectral index map shown in Fig. 16, we see 
a good correspondence between the real data and the simulation. 
Features present in the simulation also appear in the data. For in- 
stance, we see that the spectral index in the so-called Fan Region 
(i.e., near Galactic coordinates (/, b) = (90°, 20°)) is low in both 
cases, and this alone provides strong evidence for the presence 
of AME. Further, the AME spectral index is consistent with the 
spinning dust interpretation. This power-law index map may be 
used to identify particular AME regions for follow-up observa- 
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tions. Finally, we note that regions known for strong free-free 
emission, such as the Gum Nebula or Zeta Ophiuchi, have spec- 
tral indices close to -2.1 or -2.2, as expected. 

9. Conclusions 

We have tested four component separation algorithms on the 
Planck frequency maps to produce clean maps of the CMB 
anisotropics over a large area of sky. These CMB maps are used 
for studies of statistics and isotropy (Planck Collaboration XXIII 
2013), primordial non-Gaussianity (Planck Collaboration XXIV 
2013), gravitational lensing (Planck Collaboration XVII 2013), 
the ISW effect (Planck Collaboration XIX 2013), cosmic geom- 
etry and topology (Planck Collaboration XXVI 2013), searches 
for cosmic defects from primordial phase transitions (Planck 
Collaboration XXV 2013), as well as an integral part of the low- 
£ Planck likelihood (Planck Collaboration XV 2013). Two of the 
methods, one using internal foreground templates (SEVEM) and 
the other an ILC in needlet space (NILC), are non-parametric, 
extracting the CMB map by minimizing the variance of the total 
contamination. The other two methods fit models of the fore- 
grounds to clean the CMB of their emission. One fits a paramet- 
ric model in real space (C-R) and one fits a non-parametric in the 
harmonic domain (SMICA). 

All four methods have been demonstrated to work well 
both on real and simulated data, and to yield consistent results. 
Nevertheless, there are diff'erences between the methods, making 
them more or less suitable for specific applications. For instance. 
Commander -Ruler allows a joint parametric foreground estima- 
tion and CMB power spectrum estimation, with full propagation 
of foreground uncertainties to cosmological parameters, but is 
limited to a lower angular resolution than the other codes. This 
method has therefore been selected for the low-^ Planck likeli- 
hood (Planck Collaboration XV 2013) and to produce astrophys- 
ical component maps (Sect. 8), while it is sub-optimal for appli- 
cations requiring full angular resolution, e.g., gravitational lens- 
ing reconstruction or estimation of primordial non-Gaussianity. 
For these purposes, we use the three higher-resolution maps. 
We take SMICA to be the leading method, based on its supe- 
rior performance on the FFP6 simulation, where it has be shown 
to have the lowest residual foreground contamination at large 
scales and to preserve primordial non-Gaussianity. When sub- 
jecting foreground-cleaned Planck maps to scientific analysis, 
we use the other two or three maps, as appropriate, to assess the 
uncertainties inherent in the choice of methods and the assump- 
tions they make. Indeed, this is the main purpose for presenting 
four diff'erent CMB solutions to the general community. 

The CMB anisotropics are robustly recovered over a large 
fraction (73 %) of the sky and down to small angular scales, 
reaching to multipoles £ ^ 2000. We characterize the CMB maps 
with angular power spectra and cosmological parameter con- 
straints. Parameter constraints from these maps are consistent 
with those from the Planck likelihood function based on cross- 
spectra and large sky cuts (Planck Collaboration XV 2013). This 
agreement supports the robustness of both our component sepa- 
ration methodology and cosmological parameter constraints. 

The real-space parametric fits of Commander -Ruler enable 
us to characterize the diff'use Galactic foregrounds. We param- 
eterize them with a low-frequency power-law component, rep- 
resenting the sum of synchrotron, free-free, and AME emission, 
a high-frequency modified blackbody spectrum describing ther- 
mal dust emission, and a molecular CO component. Using only 
the Planck data from 30 to 353 GHz, we fit for the amplitude 
and spectral parameters of the three foregrounds and the CMB 



simultaneously at each pixel of a 40-arcmin resolution map. The 
spectral parameters are the low-frequency component power-law 
exponent and the modified blackbody emissivity power-law ex- 
ponent; the CO line ratios are spatially fixed. These parameters 
give us the source mixing matrix, which we then use in a direct 
inversion to deduce the component amplitudes at higher reso- 
lution. Through Gibbs sampling, we obtain realizations drawn 
from the full posterior distribution of possible foreground and 
CMB solutions, giving us a powerful ability to statistically char- 
acterize our results. 

Our in-depth analysis of the recovered CMB anisotropics 
is unprecedented for component separation studies, concern- 
ing both the accuracy of cosmological parameter constraints, 
and studies of early Universe physics and structure formation 
through gravitational lensing. On the other hand, the complex 
nature of the foreground emission over such a large frequency 
range limits us to the use of relatively simple methods when 
analysing Planck data alone. An extensive study in combina- 
tion with other probes of Galactic foregrounds will be presented 
in forthcoming papers. In particular, the separation of individual 
components at low frequencies requires the use of ancillary data, 
for example, from the Wilkinson Microwave Anisotropy Probe 
(WMAP) and radio surveys. 
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Appendix A: Physical parametrization 

The Commander -Ruler (C-R) approach implements Bayesian 
component separation in pixel space, fitting a parametric model 
to the data by sampling the posterior distribution for the model 
parameters. For computational reasons, the fit is performed in a 



two-step procedure: First, both foreground amplitudes and spec- 
tral parameters are found at low-resolution using Markov Chain 
Monte Carlo (MCMC)/Gibbs sampling algorithms (Jewell et al. 
2004; Wandelt et al. 2004; Eriksen et al. 2004, 2007, 2008). 
Second, the amplitudes are recalculated at high resolution by 
solving the generalized least squares system (GLSS) per pixel 
with the spectral parameters fixed to the their values from the 
low-resolution run. 

For the CMB -oriented analysis presented in this paper, we 
only use the seven lowest Planck frequencies, i.e., from 30 to 
353 GHz. We first downgrade each frequency map from its na- 
tive angular resolution to a common resolution of 40 arcmin- 
utes and repixelize at HEALPix A^side = 256. Second, we set 
the monopoles and dipoles for each frequency band as described 
in Sect. 8.1 using a method that locally records spectral indices 
(Wehus et al. 2013). We approximate the eff'ective instrumental 
noise as white with a root mean square (RMS) per pixel given 
by the Planck scanning pattern and an amplitude calibrated by 
smoothing simulations of the instrumental noise including cor- 
relations to the same resolution. For the high-resolution anal- 
ysis, the important pre-processing step is the upgrading of the 
eff'ective low-resolution mixing matrices to full Planck resolu- 
tion: this is done by repixelizing from A^side = 256 to 2048 in 
harmonic space, ensuring that potential pixelization eff'ects from 
the low-resolution map do not introduce sharp boundaries in the 
high-resolution map. 

Our model for the data, a map dy at frequency v, consists 
of a linear combination of A^c astrophysical components, plus 
instrumental noise, 

dy = ^Ft(^).A^" + ny, (A.l) 

where denotes a sky map vector containing the foreground 
amplitude map for component / normalized at a reference fre- 
quency, and F^^') is a diagonal matrix describing the spectral 
emission law for component / as a function of frequency and 
which depends on a (small) set of spectral parameters, 6. The 
CMB signal is included in the sum and, as a special case, it 
may be represented either in harmonic or pixel space, depending 
on whether the main goal of the analysis is CMB power spec- 
trum analysis or component separation. The former representa- 
tion is used for the Commander-based low-^ likelihood presented 
in Planck Collaboration XV (2013), while the latter is used for 
the foreground fits presented in Section 8 of this paper. 

Bayes theorem specifies the posterior distribution for the 
model parameters, 

P{Pl\ e\d) oc £(A\ e)P(A\ 0), (A.2) 

where £(A\ 6) is a Gaussian likelihood, and the prior P{A\ 6^) 
depends on the application. In this paper, the prior on spectral 
indices is a product of a Jeff'reys prior and physical priors, as 
detailed in Section 8.2; no priors are imposed on the foreground 
amplitudes. 

In the low-resolution Commander analysis, we exploit a 
Gibbs sampler to map out the posterior distribution (Eriksen 
et al. 2008), adopting the following minimal two-step scheme: 

A ^ P{A\e, d) (A.3) 
e ^ P(6>|A, d). (A.4) 

The first conditional distribution is a multivariate Gaussian, 
while the second distribution does not have an analytic form and 
must be mapped out numerically. 
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The high-resolution Ruler analysis maximizes the fore- 
ground amplitude conditional in Eq. A. 3 numerically by solving 
the generalized least squares system 

V fi fl 

(A.5) 

were Ny is the noise covariance matrix (assumed to be diagonal) 
of the vth channel, = F^(^) and we have neglected the differ- 
ent angular resolutions of the channels. The posterior marginal 
average for the high-resolution amplitude maps is then given by 
<A^") = / A'(e)P(e)d6 - ^ 2:^,v W(^)tDy ^ Zv WtD, a sum 

over the A^sampie samples of the spectral parameters 6. 

Once the channel weights, Wy, have been computed, pro- 
cessing a large number of simulations requires negligible com- 
putational resources. This feature has been extensively used for 
computation of the effective beam of Ruler maps: FFP6 CMB 
simulations for the 30 to 353 GHz channels are combined ac- 
cording to Wy and the effective beam transfer function is found 
as = (Q^V(C^^^w^)). Here, C^^^ is the power spectrum of the 
input simulation before convolution with the instrumental beam, 

is the HEALPix pixel window function, and the average is 
taken over the set of simulations. Missing pixels are set to 
when computing Q in the above expression. The low number 
(~ 500) of missing pixels in the data renders the impact of such 
a choice negligible at ^ < 2000. A similar procedure is used for 
defining the effective beam of the non-CMB components. 

The above algorithm produces a set of samples drawn from 
the posterior distribution, as opposed to a direct estimate of indi- 
vidual component amplitudes or spectral parameters. While this 
sample set provides a statistically complete representation of the 
posterior, it is non-trivial to visualize or to compare the distri- 
bution with external data. For convenience, we therefore sum- 
marize the distribution in terms of mean and standard deviation 
maps for each component. We emphasize, however, that the dis- 
tribution is significantly non-Gaussian, and when searching for 
features in the maps at low signal-to-noise levels, one must take 
into account the exact distribution. 

Finally, the Coimnander- Ruler confidence mask (see 
Sect. 5) is primarily defined by the product of the CG80 mask 
and the point source mask described in Sect. 4. We addition- 
ally remove any pixels excluded by the 13% Coimnander likeli- 
hood mask described by Planck Collaboration XV (2013); how- 
ever, this is almost entirely included within the CG80 mask, and 
this step therefore has very little impact on the final result. To 
complete the mask, we remove any pixels for which the high- 
resolution Ruler CMB map, smoothed to 40 arcminutes, dif- 
fers by more than 3 cr from the low-resolution Coimnander CMB 
map, which can happen due to spatial spectral variations on pixel 
scales. 

Appendix B: Internal linear combination 

NILC is a method to extract the CMB (or any component with 
known spectral behaviour) by applying the ILC technique to 
multi-channel observations in needlet space, that is, with weights 
that are allowed to vary over the sky and over the full multipole 
range. 

The ability to linearly combine input maps varying over the 
sky and over multipoles is called localization. In the needlet 
framework, harmonic localization is achieved using a set of 
bandpass filters defining a series of scales, and spatial localiza- 
tion is achieved at each scale by defining zones over the sky. 



The harmonic localization adopted here uses nine spectral bands 
covering multipoles up to ^ = 3200 (see Fig. B.l). The spatial 
localization depends on the scale. At the coarsest scale, which 
includes the multipoles of lowest degree, we use a single zone 
(no localization), while at the finest scales (which include the 
highest multipoles) the sky is partitioned into 20 zones (again, 
see Fig. B.l). 

The NILC method amounts to computing an ILC in each zone 
of each scale, allowing the ILC weights to adapt naturally to 
the varying strength of the other components as a function of 
position and multipole. A complete description of the basic NILC 
method can be found in Delabrouille et al. (2009). 

In the present work, however, there is an important differ- 
ence in the processing of the coarsest scale. Since the coarsest 
scale of the NILC filter is not localized, a plain NILC map would 
be equivalent to a pixel-based ILC for all the multipoles of that 
scale. This procedure, however, is known to be quite suscepti- 
ble to the "ILC bias" due to chance correlations between the 
CMB and foregrounds. In order to mitigate this effect, the (sin- 
gle) covariance matrix which determines the ILC coefficients at 
the coarsest scale is not computed as a pixel average, but is rather 
estimated in the harmonic domain as an average over spherical 
harmonic coefficients using a spectral weight which equalizes 
the power of the CMB modes (based on a fiducial spectrum). 
This can be shown significantly to decrease the large scale er- 
rors. 

In practice, our NILC processing depends on several imple- 
mentation choices, as follows: 

- Input channels: In this work, the NILC algorithm is applied 
to all Planck channels from 44 to 857 GHz omitting only the 
30 GHz channel. 

- Pre-processing of point sources: Identical to the SMICA pre- 
processing (see Appendix D). 

- Masking and inpainting: The NILC CMB map is actually 
produced in a three-step process. In a first step, the NILC 
weights are computed from covariance matrices evaluated 
using a Galactic mask removing about 2 % of the sky (and 
is apodized at 1°). In a second step, those NILC weights are 
applied to needlet coefficients computed over the complete 
sky (except for point source masking/subtraction), yielding 
a NILC CMB estimate over the full sky (except for the point 
source mask). In short, the weights are computed over a 
masked sky but are applied to a full sky (up to point sources). 
In a final step, the pixels masked due to point source pro- 
cessing are replaced by the values of a constrained Gaussian 
realization (inpainting). 

- Spatial localization: The boundaries of the zones used for 
spatial localisation (shown at Fig. B.l) are obtained as iso- 
level curves of a low resolution map of Galactic emission. 

- Beam control and transfer function: As in the SMICA pro- 
cessing, the input maps are internally re-beamed to a 
5 'resolution, so the resulting CMB map is automatically syn- 
thesized with an effective Gaussian beam of five arcminutes, 
according to the unbiased nature of the ILC. 

- Using SMICA recalibration: In our current implementation, 
the NILC solution uses the values determined by SMICA for 
the CMB spectrum, given in Eq. (D.6). 

Appendix C: Template fitting 

The original SEVEM algorithm produced clean CMB maps at sev- 
eral frequencies through template fitting, followed by an esti- 
mation of the CMB power spectrum from these clean maps us- 
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Table C.l: Linear coefficients, aj, and templates used to clean individual frequency maps with SEVEM. The first column lists the tem- 
plates constructed to produce clean maps. Before subtraction, the maps are smoothed to a common resolution. The (353-143) GHz 
template is constructed at the resolution of the 44 and 70 GHz channels, in order to clean the 44 and 70 GHz maps, respectively. 
For the rest of the templates, the first map used to construct the templates was filtered with the beam corresponding to the second 
map and vice versa. Note that for 100, 143 and 217 GHz channels, we give the coefficients used for the largest of the two regions 
considered in the cleaning, which covers 97 % of the sky. 



Template 


44 GHz 70 GHz 


100 GHz 


143 GHz 


217 GHz 


353 GHz 


30-70 


3.65 X 10-1 










30-44 


1.25x10-1 


-2.35 X 10-2 


2.14 X 10-2 


-1.03 X 10-1 




44-70 




1.67 X 10-1 


1.23 X 10-1 


1.76 X 10-1 




217-100 










-0.12 X lOi 


217-143 










8.99 X 10-1 


353-143 


4.05 X 10-3 9.31 X 10-3 










545-217 










9.92 X 10-2 


545-353 




5.21 X 10-3 


7.52 X 10-3 


1.84 X 10-2 




857-545 




-4.66 X 10-^ 


-6.67 X 10-5 


-1.21 X 10-4 


-5.02 X 10-4 



ing a method based on the Expectation Maximization algorithm 
(Martinez-Gonzalez et al. 2003; Leach et al. 2008; Fernandez- 
Cobos et al. 2012). From this power spectrum, a multifrequency 
Wiener-filtered CMB map was produced. For the present work, 
only the first step of the method, producing clean CMB maps 
at diff'erent frequencies, is considered. In addition, two of these 
clean maps are optimally combined to produce a final CMB map. 

The templates used for cleaning are internal, i.e., they are 
constructed from Planck data, avoiding the need for external 
data sets, which usually complicate the analyses and may in- 
troduce inconsistencies. In the cleaning process, no assumptions 
about the foregrounds or noise levels are needed, rendering the 
technique very robust. The fitting can be done in real or wavelet 
space (using a fast wavelet adapted to the HEALPix pixelization; 
Casaponsa et al. 201 1) to properly deal with incomplete sky cov- 
erage. For expediency, however, we fill in the small number of 
unobserved pixels at each channel with the mean value of their 
neighbouring pixels before applying SEVEM. 

We construct our templates by subtracting two close Planck 
frequency channel maps, after first smoothing them to a com- 
mon resolution, and converting to CMB temperature units if nec- 
essary, to ensure that the CMB signal is properly removed. A 
linear combination of the templates, tj, is then subtracted from 
(hitherto unused) map, J, to produce a clean CMB map at that 
frequency. This is done either in real or wavelet space (i.e., scale 
by scale) at each position on the sky, 

nt 

Tc(x, y) = d(x, y)-Yj ^/X-^)' (C.l) 

where nt is the number of templates. If the cleaning is performed 
in real space, the aj coeflftcients are obtained by minimizing the 
variance of the clean map, T^, outside a given mask. When work- 
ing in wavelet space, the cleaning is performed in the same way 
at each wavelet scale independently (i.e., the linear coefficients 
depend on the scale). Although we exclude very contaminated 
regions during the minimization, the subtraction is performed 
for all pixels and, therefore, the cleaned maps cover the full-sky 
(although foreground residuals are expected to be present in the 
excluded areas). 

An additional level of flexibility may also be considered: the 
linear coefficients can be fixed over the full sky, or in several 
regions. The regions are then combined in a smooth way, by 
weighting the pixels at the boundaries to reduce possible dis- 
continuities in the clean maps. 



Since the method is linear, we may easily propagate the noise 
properties to the final CMB map. Moreover, it is very fast and 
permits the generation of thousands of simulations to character- 
ize the statistical properties of the outputs, a critical need for 
many cosmological applications. The final CMB map retains the 
angular resolution of the original frequency map. 

There are several configurations possible for SEVEM, depend- 
ing the number of frequency maps to be cleaned or the num- 
ber of templates used in the fitting. Note that the production of 
clean maps at diff'erent frequencies is of great interest in order 
to test the robustness of the results. Therefore, to define the best 
strategy, one needs to find a compromise between the number 
of maps that can be cleaned independently and the number of 
templates that can be constructed. 

In particular, we have cleaned the 143 GHz and 217 GHz 
maps using four templates constructed as the diff'erence of the 
following Planck channels (smoothed to a common resolution): 
(30-44), (44-70), (545-353) and (857-545). For simplicity, the 
two maps have been cleaned in real space, since there was no 
significant improvement when using wavelets, especially at high 
latitude. In order to take into account the diff'erent spectral be- 
haviour of the foregrounds at low and high Galactic latitudes, we 
considered two independent sky regions, using diff'erent sets of 
coeflftcients (see Table C. 1 for the values of the linear coefl&cients 
for the main considered region). The first region corresponds to 
the brightest 3 % of Galactic emission, while the second region 
is defined by the remaining 97 % of the sky (see Sect. 4 for a 
detailed description of our Galactic mask construction). For the 
first region, the coefficients are actually estimated over the com- 
plete sky (we find that this is better than performing the mini- 
mization on only the brightest 3 % of the sky, where the CMB is 
very sub-dominant), while for the second region, we exclude the 
bright 3 % sky fraction, point sources detected at any frequency 
and those pixels which have not been observed in all channels. 

Note that, for consistency, we have used the same configu- 
ration (four templates, cleaning in real space, two regions) for 
the analysis of the FFP6 simulations. However, we find that this 
simple configuration produces more contaminated CMB maps 
than for the data (although the region outside the confidence 
mask still has low contamination), indicating some diff'erences 
between the foreground level in the data and in simulations. 
Therefore, conclusions derived from the FFP6 results for SEVEM 
should be taken with caution, since they are expected to provide 
overestimated residuals. 
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Multipole moment, i 




Fig. B.l: Spectral localization for NILC using nine spectral win- 
dow functions defining nine needlet scales {top panel). The 
scale-dependent spatial localization partitions the sky in one 
zone (for scale 1), two zones (for scale 2), four zones (for scale 
3), or twenty zones (for scales 5, 6, 7, 8, 9). The two-zone, 4- 
zone and 20-zone partitions are depicted in the lower panels. 



Our final CMB map was constructed by combining the 143 
and 217 GHz maps by weighting the maps in harmonic space, 
taking into account the noise level, the resolution, and a rough 
estimation of the foreground residuals of each map (obtained 
from realistic simulations). This final map has a resolution cor- 
responding to a Gaussian beam of 5' FWHM. 

Moreover, additional clean CMB maps (at frequencies 44, 
70, 100 and 353 GHz) were also produced using different com- 
binations of templates. In particular, to clean the 100 GHz map, 
we used the same templates and regions as for 143 and 217 GHz. 
This allows us to produce three (almost) independent clean 
maps, in the sense that none of the three maps to be cleaned 
is used to construct the templates. For 44, 70 and 353 GHz, dif- 
ferent combinations of templates are used and the linear coeffi- 
cients are chosen to be the same over the full sky. They are ob- 
tained by minimizing the variance of the map outside the same 
mask as that used to clean the central frequency maps on the 
largest region. The templates and the corresponding linear coef- 
ficients used for each of the considered frequencies are given in 
Table C.l. 

The SEVEM clean frequency maps have been used in analyses 
of the isotropy and statistics of the CMB (Planck Collaboration 
XXIII 2013) and to obtain cosmological constraints from the in- 
tegrated Sachs-Wolfe effect (Planck Collaboration XIX 2013). 
In particular, clean maps from 44 to 353 GHz were used for the 
stacking analysis presented in Planck Collaboration XIX (2013), 
while frequencies from 70 to 217 GHz were used for consistency 
tests in Planck Collaboration XXIII (2013). 

The confidence mask provided for SEVEM is constructed by 
combining four types of selected regions. In particular, it ex- 
cludes zones with high residuals identified through the subtrac- 
tion of SEVEM clean maps at different frequencies, as well as the 
sources detected at all Planck channels using the Mexican Hat 
Wavelet algorithm (Planck Collaboration XXVIII 2013). These 
point sources are masked with holes of varying radius, accord- 
ing to the flux of each source. The pixels that are not observed by 
all channels are also masked. Finally, to ensure that the area left 
outside the mask is statistically robust, we also exclude from the 
analysis the brightest 20 % of Galactic emission, leaving a use- 
ful area of around 76 %. This provides a conservative mask for 
CMB analysis; however, we point out that smaller masks could 
also be used in specific applications, such as the lensing potential 
reconstruction described in Sect. 7). 

Appendix D: Spectral matching 

SMICA (Spectral Matching Independent Component Analysis) 
reconstructs a CMB map as a linear combination in harmonic 
space of A^chan input frequency maps with weights that depend 
on multipole £. Given the A^chan x 1 vector x^^ of spherical har- 
monic coefficients for the input maps, it computes coefficients 
Sim for the CMB map as 

Stm^y^l-^im, (D.l) 

where the A^chan x 1 vector containing the multipole- 
dependent weights is chosen to give unit gain to the CMB with 
minimum variance. This is achieved with 



where vector a is the spectrum of the CMB evaluated at each 
channel (allowing for possible inter-channel re-calibration fac- 
tors) and is the A^chan x A/^chan spectral covariance matrix of 
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Fig. D.l: Weights, w^, given by SMICA to the input maps, after 
they are re-beamed to 5 'and expressed in ^rj, as a function of 
multipole. Top panel: linear scale; bottom panel: the absolute 
value of the weights on a logarithmic scale. 



x^yn. Taking in Eq. (D.2) to be the sample spectral covariance 
matrix of the observations. 
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(D.3) 



would implement a simple harmonic-domain ILC similar to 
(Tegmark et al. 2003). At the largest scales, we instead use a 
model, C^(^), and determine the covariance matrix to be used in 
Eq. (D.2) by fitting C^(^) to C^. This is done in the maximum 
likelihood sense for stationary Gaussian fields, that is, the best 
fit matrices, Cf(0), are obtained for 

6 = arg min ^(2^ -h 1) (QQ(6>)"^ + log det CM) • (D.4) 



Equations D.1-D.4 summarize the basic principles of SMICA; 
its actual operation depends on a choice for the spectral model 
C^(^), and on several implementation-specific details, which we 
briefly describe below. 



Fig. D.2: Contribution of each input channel to the noise in the 
SMICA map. 



We model the data as a superposition of CMB, noise and 
foregrounds. The latter are not parametrically modelled; instead, 
we represent the total foreground emission by d templates with 
arbitrary frequency spectra, angular spectra and correlations. In 
the spectral domain, this is equivalent to modelling the covari- 
ance matrices as 



(D.5) 



where Q is the angular power spectrum of the CMB, A is a 
A^chan X d matrix, is a positive definite d x d matrix, and 
is a diagonal matrix representing the noise power spectra of the 
data. The parameter vector 6 contains all or part of the quantities 
in (D.5). 

The decomposition D.5 reflects the fact that CMB, fore- 
grounds and noise are independent components of the sig- 
nal. Thus, SMICA is an ICA (independent component analysis) 
method. It operates by matching the observations C to the spec- 
tral model (D.5) using the criterion (D.4). 

The maximal flexibility in a SMICA fit of model (D.5) is ob- 
tained with all the parameters free, that is without any constraint 
on the spectrum Q, on the diagonal entries of N^, on a, or on 
A and P^. One would ideally fit all those parameters (except 
for obvious degeneracies, like that between a scale factor in a 
and the overall normalization of the CMB spectrum Q) over 
the whole multipole range. In practice, this turns out to be too 
difficult given the large dynamic range both over the sky and 
over multipoles. We resort to a pragmatic three-step approach in 
which the criterion (D.4) is minimized by first fitting a, then A, 
and finally the linear parameters Q and N^. Each fit is conducted 
over the multipole ranges and the sky fraction most appropriate 
for the parameter of interest, as follows. 

We first estimate the CMB spectral law a by fitting all model 
parameters (that is, without constraint) over a clean fraction of 
sky (fsky=^0 %) in the range 100 < ^ < 680 where the signal is 
CMB -dominated in most of the channels and the beam window 
functions are accurately known. In this fit, which is done over 
a clean part of the sky, we use a foreground emission matrix, 
A, with only four columns. From this step, we only retain the 
best fit value for vector a. In the second step, we estimate the 
foreground emissivity by fixing a to its value from the previous 
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step and fitting all the other parameters over a large fraction of 
sky (fsky=97 %) in the range 4 < ^ < 150 where the signal is 
dominated by the Galactic emission in all channels. From this 
second step, we retain the best fit value for the matrix A which, 
again, is adjusted without constraint other than having d = 6 
columns. In the last step, we fit all power spectrum parameters: 
we fix a and A to their previously found values and fit for Q and 
F{ at each £. 

Note that the first step (fitting a) amounts to re-calibrating 
the input maps on the basis of CMB anisotropics. For the maps 
in thermodynamics units, we find 

a = [0.9900,1.0000,1.0020,0.9990,1.0000, 

1 .0004, 0.9920, 1 .0457, 1 .0000] ^^'^^ 

The value at 857 GHz is not accurately recovered by SMICA, so 
we have set ags? = 1. Since the norm of a is degenerate with a 
global scale factor for the CMB angular spectrum, it can only 
be recovered by SMICA up to a scale factor. This degeneracy is 
fixed here by taking ^143 = 1. The re-calibration step could have 
been omitted since a is very close the unit vector. However, we 
found that using a improved the behavior of SMICA over using 
a = [l,...,l]. 

Before describing implementation details, we explain how 
SMICA deals with the varying resolution of the input channels, 
since the discussion thus far assumed that all input maps had 
the same resolution. Since SMICA works in the harmonic do- 
main, it is a simple matter to account for the beam transfer 
function, bi(£), of the i-th input map. The CMB sky multipoles 

contribute S£m(^ibi{i)pi(i) to the harmonic coefficient of 
the /-th map (where pi{£) is the pixel window function for the 
HEALPix map at N^^^^^- Therefore, in order to produce a final 
CMB map at 5 'resolution, close to the highest resolution of 
Planck, we only need to work with input spherical harmonics 
re-beamed to 5'; that is, to apply SMICA on vectors with en- 
tries = x^^^bsi^)/ bi(£) / pi(£), where b^i^) is a five-arcminute 
Gaussian beam function. By construction, SMICA then produces 
an CMB map with an eff'ective Gaussian beam of 5' (without the 
pixel window function). 

We now give further details on the actual implementation of 
SMICA: 

- Inputs: SMICA uses all nine Planck frequency channels from 
30 to 857 GHz, harmonically transformed up to ^ = 4000. 

- Pre-processing of point sources: SMICA is applied on input 
maps in which point sources are subtracted or masked. We 
start by fitting the PCCS point sources with SNR > 5 to 
a Gaussian shape where the source amplitude is estimated 
together with its position and a constant factor representing 
the background variance. If the fit is successful (x^ < 2), 
the fitted point source is removed from the map; otherwise 
it is masked in all channels and the hole is inpainted by a 
simple diff'usive filling process. This is done at all frequen- 
cies except 545 and 857 GHz, where all point sources with 
SNR > 7.5 are masked and inpainted. 

- Beams: When the harmonic coefficients of the input maps 
are re-beamed at 5', we do not apply exactly the expression 
•^^m ~ ■^^^^5(^)/^/(^)/A(^)i^^iitic)ned above because the fac- 
tor l/bi(£) would diverge at high £ for the lowest resolution 
input channels. That may not be a problem in infinite pre- 
cision arithmetic, but would lead to matrices C(^) with ex- 
tremely large condition numbers. Instead, we re-beam with 
the factor l//?/(^) replaced by mm(l/bi(0, 1000). The re- 
beaming of the CMB modes then is no longer perfect, but 
this is of course irrelevant because the thresholding occurs 



in a regime where the signal is completely dominated by the 
noise, so that the contribution of the corresponding channel 
is already highly attenuated by the SMICA weights (as shown 
in Fig. D.l). 

- Masking: In practice, SMICA operates on a masked sky, the 
mask being applied after the point source processing. The 
mask is obtained by thresholding a heavily smoothed version 
of the point source mask. The threshold is chosen to leave 
about 97 % of the sky. Because of the heavy smoothing, the 
mask has smooth contours and is only sensitive to large ag- 
gregates of point sources: the masked areas mostly lie in the 
Galactic plane, but include also a few bright regions like the 
Large Magellanic Cloud. 

- Inpainting: The SMICA map used in this paper has no real 
power in the masked region described above. However, for 
convenience, an inpainted SMICA map has also been pro- 
duced by replacing the masked pixels with a constrained 
Gaussian realization obtained by the method of Benoit-Levy 
et al. (2013). That map appears in Planck Collaboration I 
(2013). 

- Binning: In our implementation, we use binned spectra. 

- Processing at fine scales: Since there is little point trying to 
model the spectral co variance at high multipoles, because 
the sample estimate is sufficient, SMICA implements a simple 
harmonic ILC at ^ > 1500; that is, it appHes the filter (D.2) 
with = Q. 

- Confidence mask: A confidence mask (Fig. 3) is provided 
with SMICA, constructed in the following way. The SMICA 
CMB map is bandpass filtered through a spectral window 
v(£) = exp[-((^ - 1700)/200)V2]. The result is squared 
and smoothed at two-degree resolution, yielding a map of 
the (bandpassed) variance of the CMB map. That variance is 
corrected for the noise contribution by subtracting the vari- 
ance map for the noise obtained by the same procedure ap- 
pHed to the SMICA HRHD map. If the SMICA map contained 
only CMB and noise, the variance map would have a uni- 
form value Y.{V(0^b5(O^C(£)(2£ -\- l)/47r = 313 jjK^ over 
the sky. The confidence map is obtained by thresholding the 
noise-corrected variance map at 70yuK^. 

Viewed as a filter, SMICA can be summarized by the weights 
applied to each input map as a function of multipole. In this 
sense, SMICA is strictly equivalent to co-adding the input maps 
after convolution by specific axisymmetric kernels directly re- 
lated to the corresponding entry of w^. 

The SMICA weights used here are shown in Figure D.l (for 
input maps in units of ^rj). We see, in particular, the (expected) 
progressive attenuation of the lowest resolution channels with 
increasing multipole. Figure D.2 shows the contribution of each 
input channel to the noise in the SMICA CMB map as a function 
of multipole. The spectral noise contribution from channel / is 
simply obtained as w(£)^Ni(£), where Wi(£) is the i-th entry of 
the weight vector w(^) and Ni(£) is the angular spectrum of the 
i-th noise map. 

More details about SMICA are given in Cardoso et al. (2008), 
as well as in applications to the analysis of WMAP (Patanchon 
et al. 2005) and Archeops data (Tristram et al. 2005). An ap- 
plication to the measurement of the tensor to scalar ratio using 
CMB B-modes is discussed in Betoule et al. (2009). Within the 
Planck collaboration, SMICA is used to define the Plik high-^ 
likelihood (Planck Collaboration XV 2013), but physical mod- 
els of foreground emission are used there instead of the non- 
parametric foreground model used here. SMICA is also used to 
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Fig. E.l: Standard deviation of the four foreground-cleaned 
CMB maps from the FFP6 simulation. All maps have been 
downgraded to a HEALPix resolution of A^side = 128. 



cross-check the HFI calibration (see Planck Collaboration VIII 
(2013)). 

Appendix E: FFP6 simulation results 

We study the performance of the component separation algo- 
rithms on the FFP6 simulation described in Sect. 4, providing ad- 
ditional information beyond that in the body of the paper. Much 
of the analysis presented here mirrors that shown for the data in 
Sects. 5 and 6, allowing a direct comparison between the two 
analyses. 

E.1. CMB maps 

First, we show in Fig. E.l the standard deviation between the 
four foreground-cleaned FFP6 maps, similar to that shown in 
Fig. 5 for the data. Figure E.2 shows all pairwise differences be- 
tween the same maps, mirroring Fig. 6 for the data. These two 
plots highlight an important point concerning the FFP6 analy- 
sis already mentioned in Sect. 5, namely that in near-Galactic 
regions, where the foregrounds are important, the internal differ- 
ences between the four algorithms are larger in the FFP6 simula- 
tion than in the real data. This is due to the fact that each compo- 
nent separation algorithm has been optimized with the real data 
in terms of model definition, localization, etc. Then, the same 
models have been used for the FFP6 simulation without change. 
Only the parameters within those models are refitted to the new 
data set. This implies in fact that we expect each method to per- 
form better on the data than the simulations in terms of absolute 
residuals, to the extent that the simulation matches the real sky. 
In other words, the FFP6 simulation provides a conservative es- 
timate of the residual errors in the real data. 

E.2. Power spectrum residuals from individual components 

In Fig. E.3 we show the residual effect of some of the individual 
components on the foreground-cleaned CMB map. The thermal 
dust emission, CIB fluctuations, point sources, and noise have 
been processed individually with each algorithm. All other com- 
ponents (free-free, synchrotron, spinning dust, CO, thermal SZ, 
and kinetic SZ) are shown as a single, composite residual com- 
ponent. 



E.3. CMB power spectra and cosmological parameters 

We assess the performance of our component separation 
techniques by evaluating cosmological constraints from the 
foreground-cleaned CMB maps derived from the FFP6 simula- 
tion. 

Figure E.4 shows the estimates of the angular power spectra 
of the CMB maps. Figure E.5 compares the cosmological pa- 
rameters derived from the four foreground-cleaned CMB maps, 
CamSpec^, and Plik to the input (theoretical) parameters for 
different ^-ranges. The parameter space is defined by the same 
model applied to the real data in Sect. 6, including six ACDM 
and two foreground parameters. All deviations from input pa- 
rameters are small and within 1 cr up to ^ = 2000, verifying 
that all methods work well in this multipole range. However, for 
^max = 2500 we start to see significant shifts, e.g., for Qb^^ and 
^s- Further, the point source foreground parameter, Aps, reaches 
large values, implying that assumptions concerning the high-^ 
foreground model become important. For these reasons, we con- 
sider ^max = 2000 as the maximum recommended ^-range for 
these maps in the current data release. 

Still, the overall agreement is excellent between all codes 
and all ^-ranges. In particular, we see that differences in the band 
power spectra at high £ between the different codes are mostly 
absorbed by the two-parameter foreground model. For instance, 
the Coimnander- Ruler band power spectrum has more power at 
high £ due to noise or residual point sources, but this excess is 
well fitted by the two-parameter foreground model, and mostly 
interpreted in terms of a residual point source component; this is 
expected, given the lower angular resolution of this map. 

As mentioned above, ris and and A are to some extent sen- 
sitive to ^max- These parameters are degenerate with the fore- 
ground parameters. This may suggest that our Q foreground 
templates deviate more from the shape of the Poissonian and 
clustered component in the CMB map. This is a limitation of the 
simple foreground templates used here. To properly describe the 
foreground residuals in the reconstructed maps, we should use 
a foreground power spectrum template tailored to each method. 
For instance, such templates may be constructed by processing 
simulated foreground maps though each of the four pipelines. 
The templates are then given by the pseudo-Q of each of the 
processed foreground map. However, our analysis shows that the 
current simple model provides accurate results when restricting 
the analysis to ^max = 2000. 

Figure E.6 shows the best-fit power spectrum residuals for 
the CMB map, CamSpec and Plik relative to the input CMB 
ACDM model estimated up to ^ = 2000. These plots show that 
the residuals of the CMB map-based best-fit models are com- 
parable to the CamSpec and Plik residuals, and smaller than 
40 yuK^ for most of the £ range with larger deviations observed 
for CamSpec at ^ ~ 200. At higher £s the residuals are smaller 
than lOyuK^ for both approaches, all showing similar trends. 
Thus, both the map- and spectrum-based likelihoods recover in- 
put parameters reasonably well, with the latter yielding slightly 
larger deviations from the best-fit model of the input CMB real- 
ization. 



APC, AstroParticule et Cosmologie, Universite Paris Diderot, 
CNRS/IN2P3, CEA/lrfu, Observatoire de Paris, Sorbonne Paris 
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^ For CamSpec, ^pivot = 0.05 was adopted for this test, while all oth- 
ers, input parameters and input CMB realization included, use ^pivot = 
0.002. 
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Fig. E.2: Pairwise differences between foreground-cleaned CMB maps from the FFP6 simulation. All maps have been downgraded 
to a HEALPix resolution of A^side = 128 to show the large-scale differences. 
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Fig. E.3: Angular power spectra of residual foreground emis- 
sion in the CMB maps from the FFP6 simulation. The com- 
ponents shown are: thermal dust, cosmic infrared background 
fluctuations, point sources, CMB, noise, and the sum of all 
others. From top to bottom, the panels show the results for 
Commander-Ruler, NILC, SEVEM, and SMICA. 
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Fig. E.4: Estimates of the CMB power spectra from the 
foreground-cleaned FFP6 maps, computed by XFaster. The 
solid lines show the spectra after subtracting the best-fit model 
of residual foregrounds. The vertical dotted line shows the max- 
imum multipole {£ = 2000) used in the likelihood for fitting the 
foreground model and cosmological parameters (see Sect. E.3 
for further details). The dashed lines show the spectra before 
residual foreground subtraction. 
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Fig. E.5: Comparison of cosmological parameters derived from the FFP6 simulation using different methods. The parameters shown 
as blue, red and green points indicate results obtained with ^^ax = 1500, 2000 and 2500, respectively, and the yellow points show the 
results derived by CamSpec and Plik using cross-spectra. The black horizontal lines mark the input parameter values. The values of 
the foreground parameters are not shown for CamSpec or Plik since they use a different model. The matter power spectrum pivot 
scale was /^pivot ~ 0.002 for all likelihoods, except CamSpec for which A^pivot ~ 0.05 was used. 
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Fig. E.6: Residuals of map-based and spectrum-based best-fit 
models relative to the FFP6 simulation input ACDM spectrum, 
shown for each algorithm up to ^max = 2000. Cosmic variance is 
shown as the black dashed line. 
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